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Executive  Summary 

Accurate  tools  and  procedures  for  simulating  the  blast  effects  of  shallow  buried  explosive  devices  are 
critical  to  the  Department  of  Defense's  efforts  to  design  blast-resistant  vehicles  that  can  increase  crew 
survivability  and  counter  the  threat  from  improvised  explosive  devices.  Of  the  various  computer  codes 
available  for  modeling  explosive  effects,  the  CTH  (CHART  Squared  to  the  Three-Halves)  Eulerian 
hydrocode  from  Sandia  National  Laboratories,  NM,  is  one  of  the  most  popular  and  is  considered  a 
standard  tool  for  blast  simulations.  Several  DoD  agencies,  including  US  Army  TARDEC  is  using  CTH  as  one 
of  its  tools  for  designing  blast-resistant  vehicles.  This  project  incorporated  two  areas  of  improvement  to 
CTH  that  will  enhance  its  value  as  an  analysis  tool:  (a)  implementation  of  the  Hybrid  Elastic-Plastic  (HEP) 
Geo-material  soil  model  to  support  buried  mine  simulation,  and  (b)  procedures  for  utilizing  CTH- 
generated  blast  loads  in  LS-DYNA  structural  dynamics  simulations. 

Program  Overview 

The  code  enhancement  and  testing  phases  of  this  project  were  executed  by  Center  for  Advanced 
Vehicular  Systems  (CAVS),  Mississippi  State  University,  Starksville,  MS  (Primary  Investigator:  Dr.  Richard 
Weed).  The  R&D  effort  was  supported  by  Contract  W56HZV-08-C-0236  and  funded  as  work  directive 
WD0037  under  the  Simulation-Based  Reliability  and  Safety  (SimBRS)  program.  TARDEC  provided  the 
requirements,  test  examples  and  overall  technical  monitoring  and  guidance  to  the  code  developers  over 
the  course  of  the  project  (Technical  POC:  Dr.  Ravi  Thyagarajan) 

The  program  started  in  October  2010  and  was  completed  in  September  2012.  This  R&D  project  was 
funded  with  DoD  Ancillary  funding;  with  some  partial  funding  for  internal  TARDEC  resources  provided  by 
the  Underbody  Modeling  for  Testing  and  Evaluation  project  (Program  Manager:  Mr.  Pat  Horton, 
ARL/SLAD,  Aberdeen,  MD). 

Summary  of  Accomplishments/Results 

Both  of  the  targeted  enhancements  to  CTH  Version  9.1  described  in  the  Executive  Summary  have  been 
implemented  and  tested  in  this  project. 

Recommendations 

Both  the  HEP  geo-material  soil  model  and  the  one-way  coupler  to  LS-DYNA  have  significant  promise  in 
the  field  of  Underbody  Blast  modeling  for  lEDs.  Evaluations  should  be  continued  to  assess  and  improve 
these  features,  as  well  as  incorporate  them  into  future  versions  of  CTH,  including  the  production  version 
released  by  Sandia  National  Labs. 

These  enhancements  to  CTH  v9.1  are  available  to  any  DoD  agency  or  industry  partner  with  the  following 
restriction.  The  CTH  software  is  owned  and  maintained  by  Sandia  National  Labs,  NM  and  has  been 
designated  Export  Control  Information  (ECI)  under  the  International  Traffic  in  Arms  Regulations  (ITAR). 
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Section  1:  Introduction 

Accurate  tools  and  procedures  for  simulating  the  blast  effects  of  shallow  buried  explosive  devices  are 
critical  to  the  Department  of  Defense's  efforts  to  design  blast-resistant  vehicles  that  can  increase  crew 
survivability  and  counter  the  threat  from  improvised  explosive  devices.  Of  the  various  computer  codes 
available  for  modeling  explosive  effects,  the  CTH  (CHART  Squared  to  the  Three-Halves)  Eulerian 
hydrocode  from  Sandia  National  Laboratories  is  one  of  the  most  popular  and  is  considered  a  standard 
tool  for  blast  simulations.  TARDEC  is  using  CTH  as  one  of  its  tools  for  designing  blast-resistant  vehicles 
but  has  acknowledged  a  need  for  general  user  analysis  support  and  possible  enhancements  to  CTH  that 
will  improve  its  usability  for  TARDEC  researchers.  Two  potential  areas  of  improvement  to  CTH  that  will 
enhance  its  value  as  an  analysis  tool  for  TARDEC  are  implementation  of  material  models  for  soils  to 
support  buried  mine  simulation  and  procedures  for  utilizing  CTH  generated  blast  loads  in  LS-DYNA 
structural  dynamics  simulations. 

Although  a  wide  array  of  material  models  are  implemented  in  CTH,  it  does  not  support  the  Army 
Engineer  Research  and  Development  Center  (ERDC)  Hybrid  Elastic-Plastic  (HEP)  Geo-material  model  that 
has  proven  to  provide  accurate  simulations  of  the  response  of  a  variety  of  geo-materials  to  extreme 
loading  due  to  blast,  hypervelocity  penetration,  and  ground  shock.  Accurate  material  models  for  soils 
and  other  geo-materials  are  critical  for  analysis  of  shallow  buried  explosives  such  as  mines  and 
improvised  explosive  devices.  Since  LS-DYNA  is  the  standard  structural  dynamics  analysis  tool  at 
TARDEC,  implementing  a  procedure  for  transferring  blast  pressure  loads  generated  by  CTH  simulations 
to  LS-DYNA  provides  an  alternative  to  the  existing  LS-DYNA  Arbitrary  Lagrangian-Eulerian  (ALE)  Fluid 
Structure  Interaction  (FSI)  capability.  In  response  to  the  need  to  improve  user  support  for  CTH  and  to 
extend  the  current  capabilities  of  the  code  to  support  TARDEC  buried  explosive  simulations.  The 
Mississippi  State  University  (MSU)  Center  for  Advanced  Vehicular  Systems  (CAVS)  was  awarded  Work 
Directive  WD0037  under  the  TARDEC  Simulation  Based  Reliability  and  Safety  (SimBRS)  program. 

This  report  describes  research  and  development  activities  supported  by  Contract  W56HZV-08-C-0236 
and  funded  as  SimBRS  work  directive  WD0037  to  provide  CTH  analysis  support  and  code  enhancements. 
Subsequent  sections  of  the  report  will  describe  the  objectives  of  the  research,  the  implementation, 
verification,  and  validation  of  the  ERDC  HEP  Geo-material  model  in  CTH  version  9.1,  and  the 
development  and  implementation  of  a  coupling  procedure  that  allows  CTH  pressures  to  be  used  in  LS- 
DYNA  simulations 

1.1  Purpose  of  Project 

TARDEC  identified  the  overall  objective  of  this  Work  Directive  as  supporting  TARDEC's  efforts  to 
utilize  the  CTH  hydrocode  for  analysis  of  blast  events  on  Army  vehicles  and  to  design  efficient  systems  to 
protect  vehicle  occupants  by:  a)  providing  CTH  user  support,  b)  providing  CTH  training  as  needed  and  c) 
enhancing  CTH  usefulness  by  incorporating  TARDEC-identified  code  modifications. 


UNCLASSIFIED:  Distribution  Statement  A.  Approved  for  Public  Release 


11  |  P  a  g  e 


ENHANCEMENTS  AND  ANALYSIS  OF  CTH  SOFTWARE  FOR  UNDERBODY  BLAST 


1.2  Scope  of  the  Project  per  WD0037 

In  the  Statement  of  Work  for  this  Work  Directive,  TARDEC  identified  two  main  tasks.  A  brief  summary 
of  the  tasks  and  how  MSU  responded  to  them  follows: 


Task  1:  Provide  analysis  support  for  TARDEC  CTH  simulations 

For  Task  1,  TARDEC  requested  support  for  CTH  analysis  setup,  debugging  run  errors,  data-post 
processing  and  interpretation  of  results  from  computer  simulations.  MSU  is  committed  to  providing 
continuing  support  for  the  enhanced  version  of  CTH  developed  under  this  Work  Directive  beyond 
the  contract  completion  date. 


Task  2:  Implement  enhancements  to  CTH  to  support  blast  simulations 

For  Task  2,  TARDEC  requested  modifications  to  CTH  to  enhance  the  material  modeling  capabilities 
and  other  enhancements  with  emphasis  on  improving  CTH  analysis  capabilities  for  modeling  blast 
events  on  Army  vehicles.  Discussions  with  TARDEC  led  MSU  to  concentrate  on  two  code 
enhancements  that  TARDEC  deemed  to  be  desirable  and  of  high  priority.  The  first  code 
enhancement  requested  by  TARDEC  was  enhancing  the  material  modeling  capabilities  of  CTH  by 
incorporating  the  ERDC  HEP  Geo-Material  model  in  CTH  to  improve  the  soil  modeling  capabilities  of 
CTH  for  buried  explosive  simulations.  The  second  code  enhancement  requested  by  TARDEC  was 
implementation  of  a  one-way  "coupling"  procedure  that  would  allow  blast  pressures  generated  by 
CTH  to  be  used  in  an  LS-DYNA  structural  dynamics  simulation.  The  bulk  of  the  effort  on  this  Work 
Directive  was  devoted  to  these  two  code  enhancements. 


The  remainder  of  this  report  focuses  on  the  progress  on  the  two  code  enhancement  tasks  described 
above. 
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Section  2:  Implementation  and  verification  of  the  ERDC  Hybrid 
Elastic-Plastic  (HEP)  geo-material  model  in  CTH  9.1 

2.1  Overview  of  the  HEP  model 

The  HEP  model  was  originally  developed  by  Zimmerman  et  al.  [2]  for  modeling  ground  shock.  The 
model  was  extended  over  the  years  and  applied  to  a  variety  of  blast  and  penetration  problems  [4,  5,  6] 
ERDC  has  fit  the  model  to  a  wide  variety  of  geological  materials  ranging  from  soils  to  conventional 
strength  concretes  and  rocks.  The  model  has  many  variants  of  differing  level  of  complexity  with  the 
most  commonly  used  variants  being  the  simple  (SHEP)  and  complex  (CHEP)  models.  These  models  have 
been  validated  against  explosive  field  tests  for  a  variety  of  materials.  ERDC  maintains  a  library  of 
materials  fit  to  different  model  variants.  These  material  fits  were  generated  from  quasi-static  laboratory 
tests.  The  majority  of  the  materials  in  the  library  have  been  fit  to  the  SHEP  model.  To  date,  the  model 
has  been  implemented  in  Lagrangian  finite  element  codes  such  as  SABER,  EPIC,  and  PR0NT03D.  To  the 
best  of  our  knowledge,  the  work  described  in  this  paper  is  the  first  successful  implementation  of  the 
model  in  an  Eulerian  hydrocode. 

The  HEP  Equation  of  State  (EOS)  model  is  based  on  a  modified  Tillotson  equation  of  state  [7]  that 
calculates  an  energy  dependent  pressure-volume  response  for  a  given  material.  The  EOS  is  composed  of 
two  parts,  an  energy  independent  total  hydrostatic  response  of  constituent  solids,  air  and  water  and  an 
energy  dependent  term.  The  hydrostatic  part  of  the  HEP  EOS  consists  of  multiple  regions.  Hysteretic 
unloading  and  reloading  is  modeled  up  to  the  state  of  total  void  closure  or  "crush".  Additional  high- 
pressure  regions  and  reversible  solid-solid  phase  changes  can  be  modeled  above  void  closure.  Non- 
hysteretic  unloading  and  reloading  is  used  above  void  closure.  The  energy  independent  (hydrostatic) 
part  of  the  EOS  depends  only  on  the  compression  or  volumetric  strain  of  the  material  and  allows  for  a 
wide  variety  of  loading,  unloading,  and  reloading  conditions  in  both  compression  (taken  as  positive 
pressure)  and  tension  (negative  pressure).  The  energy  dependent  term  allows  for  liquid  phase 
transition.  For  shear  strength,  the  SHEP  model  uses  an  exponential  failure  surface  that  is  a  function  of 
pressure  and  is  defined  in  terms  of  principal  stress  differences  and  constant  Poisson  ratios  for  loading 
and  unloading  which  can  be  different  values  depending  on  the  material.  The  CHEP  model  uses  variable 
Poisson  ratios  but  a  similar  yield  surface.  For  shear  modulus,  the  SHEP  model  is  defined  by  the  current 
compression  and  local  bulk  modulus.  The  CHEP  module  uses  a  spline  fit  for  shear  modulus.  Figures  2.1a 
and  2.1b  taken  from  Reference  5  show  a  representative  pressure-compressibility  response  and  yield 
surface  for  the  simple  HEP  model.  Figure  2.1a  shows  the  various  regions  that  can  be  modeled  in  the 
hydrostatic  portion  of  the  total  equation  of  state.  A  more  detailed  description  of  the  model  can  be 
found  in  Reference  3. 
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a)  Hydrostatic  Response 


Figure  2. 1  Typical  SHEP  hydrostatic  response  and  yield  surface  (Reference  5) 


2.2  CTH  HEP  model  implementation  overview 

The  current  implementation  of  the  HEP  model  in  CTH  9.1  was  performed  in  three  phases.  In  the  first 
phase,  software  was  written  to  implement  the  Simple  HEP  EOS  and  strength  models  as  stand-alone 
routines.  These  routines  were  subjected  to  rigid  unit  testing  outside  of  CTH.  In  the  second  phase  of  the 
implementation,  the  driver  routines  and  supporting  software  required  by  CTH  were  implemented  and 
tests  were  run  using  the  CTH  Prescribed  Deformation  (PRDEF)  facility  that  is  designed  for  testing  new 
material  models  in  CTH.  The  PRDEF  tests  were  performed  for  three  HEP  soils  and  compared  with  results 
using  an  existing  CTH  soil  modeling  procedure  developed  by  Kerley  [12,  13,  14]  for  similar  soils.  In  the 
final  phase,  the  CTH  implementation  was  subjected  to  a  series  of  tests  that  compared  results  from  CTH 
simulations  with  standard  experimental  results  and  simulation  results  using  other  models  in  CTH  and 
results  from  the  EPIC  code.  The  following  sections  give  an  overview  of  the  implementation  of  the  HEP 
model  in  CTH  and  describe  results  from  the  PRDEF  tests  along  with  comparisons  with  experimental  data 
and  results  from  EPIC  code  simulations  for  three  distinct  test  cases. 

2.2.1  Implementation  of  the  EOS  routines 

The  equation  of  state  (EOS)  portion  of  ERDC  HEP  soil  model  was  implemented  using  the  Sandia  Model 
Interface  Guidelines  (MIG)  [8,9,10].  In  addition,  the  existing  Holmquist-Johnson-Cook  concrete  and 
Johnson-Holmquist  ceramics  models  [26]  were  used  as  examples.  The  CTH  implementation  began  by 
writing  new  software  based  on  the  existing  literature  [2,3]  plus  notes  provided  by  Dr.  Stephen  Akers 
(ERDC-GSL)  [11],  In  addition,  the  implementation  of  the  ERDC  HEP  model  in  PR0NT03D  was  used  as  a 
reference  in  the  development  of  the  new  routines  to  resolve  problems  found  during  testing.  However, 
none  of  the  original  ERDC  software  is  used  in  the  CTH  HEP  implementation.  Based  on  conversations 
with  Dr.  Akers,  it  was  decided  to  restrict  the  CTH  HEP  implementation  to  the  simple  or  SHEP  variant  of 
the  HEP  model  since  the  vast  majority  of  the  materials  currently  fit  to  the  HEP  model  by  ERDC  are  for 
the  SHEP  variant.  The  new  software  was  written  to  conform  to  modern  Fortran  2003  coding  standards 
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and  to  be  portable  to  a  wide  range  of  standard  conforming  compilers.  Care  was  taken  in  the 
development  of  the  core  components  of  the  model  to  insure  that  they  can  be  used  with  minor 
modifications  in  other  codes  besides  CTH. 

The  CTH  EOS  implementation  procedure  outlined  in  Hertel  and  Kerley  [9]  defined  the  CTH  EOS 
software  requirements  and  the  existing  CTH  routines  that  must  be  modified  to  implement  them.  CTH 
requires  that  two  forms  of  an  EOS  must  be  supplied  to  the  code.  Both  forms  are  required  to  return 
pressure,  sound  speed,  and  the  derivatives  of  pressure  with  respect  to  density  and  temperature  and 
energy  with  respect  to  temperature.  The  two  forms  differ  in  that  the  first  form  takes  in  density  and 
temperature  and  returns  a  specific  energy  along  with  the  other  required  values  while  the  second  form 
takes  in  density  and  specific  energy  and  returns  a  temperature.  This  requirement  introduced  three 
obstacles  to  the  CTH  HEP  EOS  implementation.  The  first  obstacle  is  that  the  ERDC  HEP  EOS  is  cast  in  a 
form  that  is  dependent  on  internal  energy.  Therefore,  a  procedure  to  extract  a  temperature  that  is 
consistent  with  CTH  is  required.  The  most  obvious  way  is  to  assume  (probably  erroneously)  a  change  in 
internal  energy  is  a  linear  function  of  a  change  in  temperature  and  the  two  are  related  by  a  specific  heat 
at  constant  volume.  This  is  the  approach  taken  by  both  the  Holmquist-Johnson-Cook  (HJC)  concrete  and 
Johnson-Holmquist  ceramic  models  in  CTH  and,  until  a  more  appropriate  relationship  is  found,  is  the  one 
used  for  this  initial  implementation  of  CTH  HEP.  However,  this  choice  leads  to  the  second  obstacle, 
specifically;  how  to  define  representative  values  for  the  specific  heat  of  soil.  This  information  is  not 
included  in  the  ERDC  HEP  material  library.  The  approach  used  to  define  an  appropriate  specific  heat  is 
given  in  the  following  section.  The  third  obstacle  was  the  difference  in  the  reference  conditions  for 
thermodynamic  variables  used  by  HEP  and  CTH.  The  HEP  EOS  model  is  referenced  to  gage  pressure  but 
CTH  uses  an  absolute  pressure.  Since  the  HEP  model  contains  an  energy  dependent  term,  an  "energy 
shift"  has  to  be  applied  to  the  absolute  energy  value  computed  by  CTH  to  insure  that  CTH  will  return 
consistent  values  of  initial  pressure  and  temperature.  The  energy  shift  is  computed  using  the  specific 
heat  and  a  value  for  room  temperature  (298  degrees  K).  A  similar  procedure  for  defining  a  shifted 
reference  energy  state  is  used  in  the  existing  CTH  Johnson-Holmquist  ceramic  model  to  overcome 
problems  associated  with  very  low  or  even  negative  energies. 

Once  the  temperature  -  energy  relationship  was  settled  on  and  a  procedure  implemented  for 
defining  a  specific  heat,  two  driver  routines  were  written  that  take  in  the  required  dependent  variables 
and  then  call  the  base  HEPEOS  routine  to  return  pressure  and  the  other  required  values.  Since  the  HEP 
models  use  Mbar  as  the  unit  of  pressure,  the  pressure  and  energy  terms  are  scaled  on  entry  to  and  exit 
from  the  driver  routines  to  take  in  and  return  the  default  CTH  force  units  of  dynes.  The  next  step  was  to 
define  the  required  extra  or  history  variables  that  must  be  carried  over  from  the  previous  time  step.  The 
SHEP  EOS  implementation  only  requires  one  history  variable,  the  maximum  value  of  the  compression  or 
volumetric  strain  for  each  cell.  The  strength  model  requires  at  least  one  more  extra  variable,  the 
derivative  of  pressure  with  respect  to  compression.  For  the  initial  implementation  we  are  only  defining 
the  maximum  compression  over  time.  Next  a  material  number  based  on  the  range  of  available  material 
numbers  in  CTH  was  selected.  For  the  current  CTH  HEP  implementation,  the  both  the  Equation  of  State 
and  Elastic-Plastic  model  material  numbers  was  set  to  27.  The  last  tasks  before  final  implementation  of 
the  model  were  to  create  a  data  check  routine  to  initialize  computed  constants  and  check  input  data. 
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write  an  extra  variable  routine  that  is  used  to  register  the  extra  variable  with  the  CTH  database  and 
modify  a  relatively  small  group  of  CTH  routines  to  include  the  new  model.  The  only  other  modifications 
were  the  addition  of  a  module  routine  that  holds  the  global  variables  for  the  HEP  model  along  with 
some  support  routines  and  updating  the  existing  CTH  make  files  to  include  the  new  software  in  the 
compilation.  A  more  detailed  description  of  the  new  subroutines  that  make  up  the  HEP  EOS  model  and 
the  changes  made  to  existing  CTH  routines  is  given  in  Appendix  A. 

2.2.2  Strength  model  implementation 

The  constitutive  model  was  implemented  in  CTH  as  a  "traditional"  total  stress  model  that  is  based  on 
Von  Mises  J2  plasticity.  The  HEP  routines  supply  CTH  with  a  value  of  the  material  yield  surface  based  on 
the  deviatoric  principle  stresses  and  an  estimate  of  shear  modulus.  CTH  weights  these  values  with  the 
volume  fraction  of  the  material  and  sums  these  values  with  the  values  for  other  materials  in  a  cell  that 
can  support  material  strength  to  produce  a  cell  value  for  shear  modulus  and  yield  strength.  Following 
the  Johnson-Holmquist  ceramic  model  implementation,  a  set  of  seven  include  files  were  developed  that 
define  the  number  of  variables  used  in  the  model,  the  extra  variables  required  by  the  model,  and  a 
driver  routine  that  calls  the  HEP  yield  surface  and  shear  module  routines.  As  with  the  EOS  component  of 
the  HEP  model,  a  more  detailed  description  of  the  code  and  changes  to  existing  CTH  software  is  given  in 
Appendix  A.  Although  soils  in  general  have  little  strength  in  tension,  the  constitutive  model  plays  an 
important  part  in  the  values  of  impulse  and  crater  geometry  generated  by  the  HEP  model  for  soils  with 
low  adhesion  such  as  dry  sand. 

2.2.3  The  GenCTHInput  utility 

In  addition  to  the  CTH  software,  an  utility  (named  GenCTHInput)  was  written  to  translate  HEP  model 
inputs  provided  by  ERDC  in  an  input  file  format  used  in  EPIC  and  PR0NT03D  into  a  format  that  can  be 
read  by  the  standard  CTH  input  routines.  GenCTHInput  is  a  standalone  code  written  in  Fortran  2003  and 
will  be  provided  with  the  enhanced  version  of  CTH  9.1.  This  utility  was  written  to  allow  for  the  addition 
of  new  HEP  materials  not  included  in  the  current  ERDC  HEP  library.  In  addition  to  translating  the  EPIC 
and  PR0NT03D  format  files  into  a  format  the  CTH  can  read,  GenCTHInput  computes  the  value  of 
material  specific  heat  using  the  procedure  described  in  the  next  section.  The  GenCTHInput  utility  will  be 
described  in  more  detail  in  the  Users  Manual  that  will  accompany  release  of  the  enhanced  version  of 
CTH. 

2.2.4  Definition  of  soil  specific  heat  values 

A  survey  of  existing  literature  on  the  thermal  properties  of  soils  revealed  that  the  most  common 
procedure  for  estimating  the  volumetric  heat  capacities  or  specific  heats  of  soils  was  introduced  by  De 
Vries  [16].  De  Vries  computes  the  volumetric  heat  capacity  of  a  soil  as  a  sum  of  the  heat  capacities  of  the 
individual  soil  constituents,  the  heat  capacity  of  water  scaled  by  the  gravimetric  water  content,  the  heat 
capacity  of  air  based  on  the  void  fraction  of  air  in  the  soil,  and  the  heat  capacity  of  the  organic  matter  in 
the  soil.  Several  researchers  [17,18,19]  use  a  reduced  form  of  this  approximation  that  involves  a  mean 
value  for  the  dry  soil  constituents  (usually  taken  to  be  in  the  range  of  0.73  kJ/kg/K  to  0.837  kJ/kg/K),  the 
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heat  capacity  of  water  (4.180  kJ/km/K),  the  dry  bulk  density  of  the  solids  and  the  soil  water  content  in 
terms  of  weight  fraction.  The  form  of  the  equation  for  the  heat  capacity  then  becomes: 

pCv=pb(cs+CwW)  (Eq.  1) 

where  pcv  is  the  volumetric  heat  capacity  of  the  soil,  pb  is  the  dry  bulk  density  for  the  soil,  cs  is  the 
specific  heat  of  the  dry  soil  solids  (taken  to  be  0.837  kJ/kg/K  for  this  report),  cw  is  the  specific  heat  of 
water  (4.184  kJ/kg/K  where  the  density  of  water  is  1  gm/cc)  and  w  is  some  approximation  of  the 
gravimetric  moisture  content  of  the  soil  (kg/kg).  This  can  be  either  the  water  content  in  percent 
measured  by  experiment  or  the  ratio  of  the  dry  to  wet  material  densities.  The  option  to  use  this 
equation  or  input  a  specific  heat  value  directly  is  included  in  the  GenHEPInput  utility  which  will  convert 
the  results  into  the  CTH  units  for  specific  heat  of  ergs/g/ev  where  ev  is  an  electron  volt  which  is 
equivalent  to  11604.5  degrees  Kelvin.  For  the  soils  shown  in  Table  2.1,  this  procedure  leads  to  values  of 
Cv  on  the  order  of  1.1E11  to  1.7E11  ergs/g/ev.  These  values  are  consistent  with  Cv  values  of  other 
materials  in  the  CTH  material  database.  An  important  point  to  note  is  that  the  specific  heat  value 
computed  by  this  procedure  is  primarily  a  function  of  the  moisture  content  of  the  soil.  Typical  values  for 
the  wet  and  dry  densities  along  with  water  content  and  percent  air  voids  are  given  in  Table  2.1  for  the 
six  HEP  materials  used  in  this  report. 


Material 

Description 

ERDC  Material 

Classification 

Soil 

Classification 

Air  Voids, 

% 

Wet 

Density, 

Mg/MA3 

Dry 

Density, 

Mg/MA3 

Water 

Content,  % 

Gravelly  Clay 

WCLAY1 

CL 

3.9 

1.96 

1.59 

23.7 

Clayey  Sand 

ISOIL1 

SC 

12.0 

2.02 

1.80 

12.1 

Concrete 

Sand 

DSOIL1 

SP 

26.5 

1.86 

1.79 

4.0 

Clay 

IMD  WET  CLAY 

CL 

3.4 

1.97 

1.59 

24.0 

Silty  Sand 

IMD  INT SOIL 

SC 

10.9 

2.01 

1.84 

12.7 

Dry  Sand 

IMD  DRY  SAND 

SP 

29.8 

1.77 

1.70 

3.8 

Table  2.1.  Properties  of  Six  ERDC  HEP  soils 
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2.3  Verification  testing  of  the  HEP  model  implementation 

2.3.1  Unit  testing  of  the  CTH  HEP  model  software 

A  unit  test  driver  for  the  new  HEP  software  was  written  to  exercise  the  software  outside  of  CTH.  This 
testing  was  a  critical  part  of  the  development  process  because  it  helped  identify  bugs  that  would  have 
been  hard  to  diagnose  in  a  standard  CTH  run.  The  driver  allowed  the  EOS,  yield  surface  and  shear 
modulus  components  to  be  tested  using  either  a  user  defined  strain  history  or  input  from  a  CTH  run.  The 
user  defined  strain  option  allowed  for  the  definition  of  multiple  unloading-reloading  paths.  The 
objectives  of  the  tests  were  to  reproduce  Figures  5  and  6  from  Reference  4  and  compare  results  with 
the  output  from  the  ERDC  HEP  routines.  These  results  are  shown  in  Figures  2.2a  and  2.2b  for  the 
gravelly  clay,  clayey  sand,  and  concrete  sand  data  HEP  materials.  The  exact  correlation  between  the 
new  routines  and  the  old  ERDC  routines  indicated  that  the  initial  implementation  of  the  hydrostatic  part 
of  the  EOS  and  the  yield  surface  was  correct. 


Volumetric  Strain  (percent)  Mean  Normal  Stress  (MPa) 


a)  Hydrostatic  Response  b)  Yield  Surface 

Figure  2.  2:  Comparison  of  CTH  HEP  and  ERDC  software  results  for  hydrostatic  response  and  yield  surface 


2.3.2  CTH  Prescribed  Deformation  (PRDEF)  test:  Results  and  Discussion 

A  next  set  of  tests  were  run  following  a  suggestion  by  Dr.  Akers  to  use  the  CTH  Prescribed 
Deformation  (PRDEF)  option  to  mimic  a  uniaxial  compression  test  in  CTH.  The  PRDEF  option  allows  the 
user  to  run  CTH  with  a  user  specified  strain  time  history  using  different  strain  measures.  PRDEF  is 
designed  for  testing  CTH  EOS  and  strength  models  separate  from  the  other  parts  of  the  code.  The  strain 
history  used  for  these  tests  was  taken  from  a  CTH  example  problem  and  is  shown  in  Figure  2.3.  Tests 
were  run  using  the  HEP  model  fits  for  the  WCLAY1  and  DSOIL1  materials  and  compared  with  results  for 
similar  soils  using  inputs  developed  by  Kerley  [12.  13.  14].  The  Kerley  reports  outline  a  procedure  for 
using  existing  CTH  SESAME  EOS  tables  [8]  for  soils  with  various  water  contents  generated  by  Kerley  to 
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model  some  of  the  soils  in  the  HEP  database.  These  tests  resulted  in  the  pressure  versus  strain  and 
pressure  versus  time  results  shown  in  Figures  2.4a-2.4d  for  the  WCLAY1  and  DSOIL1  materials.  The  CTH 
HEP  implementation  is  seen  to  agree  well  with  the  Kerley  model  for  the  wetter  material  (WCLAY1)  but 
to  yield  substantially  different  results  for  the  dry  sand  material  (DSOIL1).  However,  it  has  been 
demonstrated  by  comparisons  with  experimental  results  that  the  Kerley  fit  is  not  a  good  model  for  this 
material.  These  results  demonstrated  that  the  CTH  HEP  implementation  delivered  equal  or  better 
material  response  for  an  arbitrary  strain  history  when  compared  to  the  Kerley  models. 

CTH  Prescribed  Deformation  Test  Strains 


Uniaxial  Strain  vs  Time 


Figure  2.  3:  Strain  vs  time  history  for  PRDEF  tests 


Pressure-Strain  Response  for  cth  Prescribed  Deformation  Test  Pressure  vs  Time  for  CTH  Prescribed  Deformation  Test 

Comparison  of  CTH  HEP  and  Kerley  SESAME  table  results  „„  Comparison  of  CTH  HEP  and  Kerley  SESAME  table  results 


2.4a)  WCLAY1  P  vs  Strain 


2.4b)  WCLAY1  P  vs  Time 
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Pressure-Strain  Response  for  CTH  Prescribed  Deformation  Test 
Comparison  of  CTH  HEP  and  Kerley  Sesame  Table  Results 


Pressure  vs  Time  for  CTH  Prescribed  Deformation  Test 

Comparison  of  CTH  HEP  and  Kerley  Sesame  Table  Results 


2.4c)  DS0IL1  P  vs  Strain  2.4d)  DS0IL1  P  vs  Time 

Figure  2. 4:  PRDEF  Hydrostatic  Response  of  CTH  HEP  and  Kerley  Soil  Models  for  WCLAY1  and  DS0IL1  materials 


2.3.3  Pioneer  buried  explosive  simulations 

Initial  pioneer  testing  with  the  CTH  HEP  implementation  was  performed  using  the  test  problem 
described  in  Kerley's  reports  along  with  a  similar  test  case  used  by  Littlefield  [15]  to  compare  results 
from  the  EPIC  HEP  implementation  with  CTH  runs  using  Kerley's  model.  Both  test  cases  simulate  blast 
effects  on  a  witness  plate  suspended  over  a  shallow-buried  explosive  charge.  Cases  for  three  HEP  soils 
were  run  with  the  Kerley  inputs  as  well  as  the  new  HEP  model  and  the  results  were  compared.  The  CTH 
HEP  results  were  also  compared  with  results  presented  by  Littlefield  that  compares  results  using  the 
Kerley  SESAME  model  in  CTH  results  with  the  EPIC  code  using  the  HEP  model  for  a  witness  plate 
simulation  similar  to  the  Kerley  tests. 

2.3. 3.1  Kerley  test  case:  Results  and  Discussion 

The  initial  tests  for  the  CTH  implementation  were  performed  using  the  test  case  described  by  Kerley 
[13]  for  the  three  soils  mentioned  previously.  The  test  case  consists  of  a  20  cm  thick  witness  plate  of 
RHA  steel  with  a  diameter  of  276  cm  suspended  46  cm  above  a  soil  mass  that  is  124  cm  deep  and  368 
cm  wide.  A  9217  g  (20.32  lb)  charge  of  TNT  is  buried  in  the  soil  at  a  depth  of  burial  of  16  cm.  The  TNT  is 
taken  to  be  a  cylinder  with  a  radius  of  15  cm  and  a  height  of  8  cm.  The  TNT  is  modeled  in  these 
simulations  using  the  JWL  equation  of  state.  Kerley  suggests  using  a  procedure  that  allows  the  SESAME 
table  values  for  the  TNT  detonation  products  to  be  used  in  place  of  JWL  in  order  to  avoid  potential 
problems  with  JWL  at  low  pressures  but  JWL  was  deemed  sufficient  for  these  tests.  The  RHA  steel  plate 
is  modeled  using  the  SESAME  EOS  for  RHA  steel  along  with  the  CTH  elastic-plastic  model  and  the 
Johnson-Cook  fracture  model  with  a  strength  in  tension  of  -3.8E10  dynes/cm2.  All  tests  were  run  in  2D 
cylindrical  coordinates  using  a  constant  mesh  spacing  of  1  cm  in  both  directions.  Half-symmetry  was 
assumed  to  reduce  the  size  of  the  problem  space. 

The  three  soils  were  modeled  using  the  SESAME  table  inputs  generated  by  Kerley  [13].  Kerley's  model 
takes  SESAME  tables  that  were  generated  for  soils  with  various  water  contents  and  then  adjusts  them 
with  the  CTH  p-alpha  porosity  model  to  match  the  dry  density,  porosity,  and  water  content  of  the  three 
soils  used  in  these  tests.  The  model  attempts  to  match  as  closely  as  possible  Hugoniot  data  obtained 
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from  ERDC  and  other  sources.  Since  the  HEP  strength  model  was  not  implemented  at  the  time  of  these 
initial  tests,  the  elastic-plastic  and  CTH  geo  model  inputs  suggested  by  Kerley  were  used  to  model  the 
shear  response  of  the  three  soils  for  both  the  Kerley  input  and  the  CTH  HEP  simulations. 

Simulations  for  all  three  soils  were  performed  with  simulation  times  that  ranged  from  3  ms  to  5  ms. 
The  3ms  simulations  were  used  to  compare  with  the  material  plots  given  by  Kerley  for  the  three  cases. 
The  5  ms  simulations  were  used  to  compare  impulse  time  histories  for  RHA  plate.  Impulse  for  these 
correlations  was  taken  to  be  equivalent  to  the  total  Y  momentum  of  the  RHA  material.  The  three  cases 
using  the  Kerley  SESAME  input  were  run  first  followed  by  the  same  cases  using  the  CTH  HEP  model.  The 
only  change  in  the  inputs  was  to  replace  the  SESAME  EOS  input  for  each  soil  with  its  equivalent  CTH  HEP 
input.  Figures  2.5a-2.5f  compare  the  combined  material  and  pressure  plots  for  the  Kerley  model  and 
CTH  HEP  at  3.0  milliseconds  for  each  of  the  three  soils.  The  material  plots  are  the  left  half-plane  of  the 
figures.  Comparison  of  the  Kerley  results  with  the  results  given  in  his  report  and  with  the  CTH  HEP 
results  revealed  some  interesting  differences.  First,  the  Kerley  results  for  CTH  V9.1  do  not  match  the 
results  for  the  older  CTH99  version  used  in  his  report.  In  particular,  there  appears  to  be  more  soil 
material  in  regions 


2.5a.)  Kerley  -  Clayey  Sand 


2.5b.)  HEP -Clayey  Sand 
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2.5c)  Kerley  -  Gravelly  Clay 


2.5d.)  HEP -Gravelly  Clay 
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2.5e.)  Kerley  -  Concrete  Sand  2.5f.)  HEP  -  Concrete  Sand 

Figure  2.  5:  Comparison  of  Kerley  and  HEP  crater  results  for  Kerley  test  case  at  3.0  ms 


where  the  older  version  of  CTH  shows  detonation  products.  However,  the  CTH  HEP  results  appear  to 
more  closely  match  the  CTH99  results  than  the  results  using  CTH  V9.1  with  Kerley's  inputs.  Littlefield 
[15]  saw  similar  differences  in  the  simulations  he  ran  to  compare  with  the  EPIC  Lagrangian  code.  A  set  of 
CTH  material  discards  that  removed  some  soil  material  at  very  low  pressures  were  required  to  yield 
material  plots  that  matched  the  EPIC  results.  It  is  felt  we  are  encountering  a  similar  problem  with  the 
Kerley  test  case.  The  effect  of  removing  the  low-pressure  material  will  be  illustrated  in  a  following 
section  where  a  comparison  with  results  generated  by  one  of  Littlefield's  test  cases  will  be  presented. 

As  observed  by  Kerley,  the  concrete  sand  material  which  has  the  lowest  water  content  of  the  three 
soils  (see  Table  2.1)  and  the  highest  yield  strength  has  the  lowest  amount  of  lateral  containment  of  the 
detonation  products  by  the  soil.  The  moister  soils  (designated  gravelly  clay  and  clayey  sand)  show 
considerable  lateral  containment  of  the  detonation  products.  The  effect  of  the  increased  amount  of 
lateral  containment  can  be  seen  in  the  comparison  of  impulse  versus  time  for  the  three  soils  given  in 
Figure  2.6.  The  moister  soils  are  seen  to  generate  much  higher  impulse  on  the  witness  plate.  The  excess 
amount  of  soil  material  in  the  CTH  V9.1  results  makes  estimates  and  comparison  of  crater  dimensions 
difficult.  However,  a  visual  comparison  of  the  CTH  HEP  results  with  the  CTH99  results  in  Kerley's  report 
shows  that  CTH  HEP  is  giving  similar  crater  sizes.  In  general,  the  correlations  for  the  gravelly  clay  and 
clayey  sand  look  very  good  with  the  clayey  sand  case  giving  an  excellent  correlation  for  impulse.  The 
concrete  sand  material  correlations  can  be  characterized  as  fair  to  good  at  best. 
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Correlation  of  CTH  HEP  and  Kerley  Soil  Models  for  Three  Soils 

Impulse  vs  Time  for  RHA  Witness  Plate 


Figure  2.  6:  Correlation  of  impulse  on  a  witness  plate  using  CTH  HEP  and  Kerley  soil  models  for  three  test  soils 


Another  interesting  result  of  these  tests  is  that  the  CTH  HEP  EOS  appeared  to  be  more  stable  than 
Kerley's  model.  The  clayey  sand  case  would  not  run  with  Kerley's  inputs  without  using  material  discards 
to  remove  material  with  negative  energies  and  densities.  All  the  cases  run  with  the  CTH  HEP  model  ran 
to  the  specified  stop  time  without  discards.  However,  subsequent  tests  did  require  discards  for  stable 
simulations  with  the  HEP  model. 

2.3. 3.2  Littlefield  test  case:  Results  and  Discussion 

Dr.  David  Littlefield  (UAB)  compared  of  Kerley's  model  for  the  three  soils  used  in  this  report  with 
results  from  the  EPIC  code  [15].  Dr.  Littlefield  provided  input  decks  for  his  Kerley  model  comparisons 
along  with  impulse  data  from  the  EPIC  code.  Littlefield's  test  case  was  run  using  the  WCLAY1  and 
DSOIL1  HEP  inputs  and  compared  with  his  results  for  Kerley's  model  and  EPIC.  The  material  discards 
criterions  provided  in  Dr.  Littlefield's  input  were  used  to  provide  a  better  visual  estimate  of  the  crater 
depth  and  diameter  generated  by  the  CTH  HEP  model  and  to  provide  data  to  compare  with  the  EPIC  and 
CTH  results  using  Kerley's  input.  Dr.  Littlefield's  test  case  differs  slightly  from  the  Kerley  test  case, 
primarily  in  the  extent  of  the  computational  domain  and  the  depth  of  burial  of  the  charge.  For  these 
tests  the  charge  is  buried  10  cm  below  the  soil  surface  instead  of  16  cm.  Also,  the  Johnson-Cook  model 
for  steel  was  used  for  the  strength  model.  A  tensile  fracture  stress  of  -3.018E9  dynes/cm2  was  used  for 
the  CTH  fracture  model.  Results  were  generated  using  the  CTH  HEP  model  with  the  material  strength 
model  used  by  Kerley  and  the  HEP  strength  model.  The  initial  material  configuration  for  Littlefield's  test 
case  is  shown  in  Figure  2.7.  Material  response  and  pressure  at  3  ms  are  shown  in  Figures  2.8a  and  2.8b 
for  the  WCLAY1  material.  The  computed  impulse  for  the  two  models  is  compared  with  EPIC  results  in 
Figure  2.9.  The  results  in  Figure  2.9  illustrate  both  the  minimal  effect  of  the  CTH  strength  model  for  the 
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WCLAY1  material  and  the  improved  accuracy  of  the  Eulerian  solution  over  the  Lagrangian  solution 
generated  by  EPIC  when  the  material  behaves  more  like  a  fluid  than  a  solid.  These  results  are  contrasted 
with  the  results  shown  in  Figures  2.10  and  2.11  for  the  dryer  DSOIL1  material.  The  computed  impulse 
from  the  CTH  HEP  model  and  the  results  from  EPIC  are  in  excellent  agreement.  The  Kerley  inputs  for  this 
case  under-predicts  the  impulse.  With  material  discards  in  place,  the  resulting  material  plots  for  both 
materials  show  good  agreement  for  the  depth  and  diameter  of  the  crater. 


Figure  2.  7:  Initial  material  configuration  and  pressure  for  Littlefield  test  case 


2.8a)  CTH  HEP 


2.8b)  Kerley  Model 


Figure  2.  8:  Littlefield  test  WCLAY1  material  response  and  pressure  at  3ms  for  CTH  HEP  and  Kerley  models 
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Impulse  vs  Time  for  WCLAY1  (gravelly  clay) 


Comparison  of  CTH  HEP,  Kerley  model  and  EPIC  HEP  for  Littlefield  case  results 


Figure  2.  9:  Correlation  of  Littlefield  test  WCLAY1  impulse  for  CTH  HEP  and  Kerley  models  with  EPIC  results 


2.10a)  CTH  HEP 


2.10b)  Kerley  Model 


Figure  2. 10:  Littlefield  test  DSOIL1  material  response  and  pressure  at  3ms  for  CTH  HEP  and  Kerley  models 
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Impulse  vs  Time  for  DSOIL1  (concrete  sand) 


Comparison  of  CTH  HEP,  Kerley  model,  and  EPIC  HEP  Littlefield  case  results 


Figure  2. 11:  Correlation  of  Littlefield  test  DSOIL1  impulse  for  CTH  HEP  and  Kerley  models  with  EPIC  results 


2.4  ERDC  Impulse  Measurement  Device  (IMD)  simulations 

The  final  sets  of  validation  tests  for  the  HEP  model  were  run  using  the  three  soils  tested  in  the  ERDC 
Impulse  Measurement  Device  (IMD)  [20].  The  majority  of  the  results  generated  for  these  tests  used  a 
simplified  version  of  the  IMD  geometry  in  which  the  mass  of  the  piston  and  guide  rails  of  the  IMD  device 
were  modeled  by  increasing  the  height  of  the  target  plate  component  to  match  the  overall  mass  of  the 
system.  The  simplified  geometry  was  adopted  to  limit  the  extent  of  the  computational  domain  so  that 
runs  could  be  made  on  relatively  course  grids  to  improve  turnaround  times  for  the  tests.  Simulation  of 
the  explosive  process  was  limited  to  the  CTH  programmed  burn  model  and  the  JWL  equation  of  state.  In 
addition,  only  the  soil  in  the  test  bed  was  modeled.  A  more  detailed  CTH  model  of  the  IMD  experiment 
was  obtained  from  Dr.  Greg  Bessette  (ERDC/GSL)  [25]  late  in  the  program  that  includes  part  of  the 
piston  assembly  and  models  both  the  test  bed  soils  and  the  surrounding  in  situ  soils.  The  CTH  input  file 
provided  by  Dr.  Bessette  contained  a  fit  of  the  IMD  WET  CLAY  material  obtained  using  the  Kerley 
SESAME  tables  and  P-alpha  porosity  corrections.  The  explosive  detonation  process  is  modeled  using  a 
reactive  burn  model  instead  of  JWL.  Due  to  time  and  budget  constraints  only  the  WET  CLAY  case  was 
run  using  Bessette's  detailed  model  of  the  IMD  test  configuration. 

For  all  IMD  test  cases,  the  computed  impulse  and  crater  dimensions  for  the  buried  charge  cases  were 
compared  with  the  experimental  values  and  the  EPIC  calculations  performed  by  Moral  et  al.  The  HEP 
material  model  data  for  the  three  soils  and  the  EPIC  impulse  results  were  obtained  from  ERDC.  The 
following  sections  present  an  overview  of  the  IMD  experiments,  the  three  soils  used  in  the  test  and 
results  obtained  from  the  simulations. 


UNCLASSIFIED:  Distribution  Statement  A.  Approved  for  Public  Release 


26  |  P  a  g  e 


ENHANCEMENTS  AND  ANALYSIS  OF  CTH  SOFTWARE  FOR  UNDERBODY  BLAST 


2.4.1  Overview  of  the  ERDC  IMD  facility 

The  ERDC  IMD  facility  shown  in  Figure  2.12  was  designed  and  built  by  ERDC  under  TARDEC  funding  to 
obtain  data  on  the  aboveground  environment  created  by  the  detonation  of  shallow-buried  explosives. 
The  facility  was  designed  to  obtain  data  on  overpressure,  impulse,  and  ground  shock  for  different  soil 
types.  A  complete  description  of  the  IMD  facility,  the  testing  procedures  and  instrumentation,  and  the 
test  conducted  with  the  facility  is  given  in  Reference  20.  ERDC  developed  three  test  soils  for  the  IMD 
experiments  that  were  based  on  the  WCLAY1,  ISOIL1,  and  DSOIL1  materials  [21,22,23].  The  IMD 
experiments  included  above  ground,  partially  buried  (surface  tangent),  and  shallow  buried  (Depth  of 
Burial  of  4  inches)  detonations  in  three  soils  of  different  type  and  water  content.  In  addition  to  impulse 
measurements,  side-on  pressure  and  ground  shock  data  was  obtained.  The  results  presented  in  the 
following  sections  will  concentrate  on  the  shallow  buried  cases  since  they  were  deemed  to  be  of  the 
most  immediate  interest  to  TARDEC.  The  ERDC  IMD  experiments  represent  the  most  complete  set  of 
data  on  the  effect  of  soil  types  and  explosive  Depth  of  Burial  (DOB)  to  date  and  provide  an  invaluable 
resource  for  validation  of  soil  models  used  in  blast  simulations.  The  distance  from  the  bottom  of  the 
plates  to  the  top  of  the  charge  was  maintained  at  20  inches  for  all  the  test  configurations. 


Figure  2.  12:  Schematic  of  ERDC  Impulse  Measurement  Device  (Reference  20) 


2.4.2  Description  of  the  three  IMD  soils 

The  three  soils  tested  in  the  IMD  facility  consisted  of  a  wet  clay  material,  and  intermediate  soil 
described  as  "silty  sand",  and  a  dry  sand  material.  Table  2.1  presents  the  material  densities  and  water 
content  of  the  three  soils.  ERDC  characterized  these  three  materials  and  fit  the  material  properties  to 
the  Simple  HEP  (SHEP)  model.  These  material  fits  were  used  by  Moral  et  al.  [5]  in  EPIC  simulations  that 
compared  the  computed  impulse  with  the  experimental  values.  These  fits  were  obtained  from  Dr.  Moral 
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and  used  in  the  CTH  HEP  simulations  without  modification.  CTH  input  files  for  the  three  soils  were 
generated  using  the  GenHepinput  tool.  As  shown  in  the  results  presented  in  the  following  sections,  the 
soil  type  (in  particular  the  water  content)  along  with  the  DOB  has  a  dramatic  effect  on  the  amount  of 
impulse. 

2.4.3  Simulation  of  the  ERDC  IMD  experiments  with  a  simplified  geometry 

A  set  of  simulations  were  run  using  the  CTH  HEP  model  with  the  goal  of  comparing  the  results  for 
impulse  and  crater  geometry  for  a  shallow  buried  explosive  with  Ehrgot's  experimental  results  and 
Moral  et  al.'s  EPIC  simulations.  Simulations  were  run  using  a  5  Ibm  charge  of  C4  at  a  DOB  of  4  inches  for 
the  three  soils  tested  in  the  IMD,  as  reported  in  Reference  20.  The  CTH  simulations  used  a  2D  cylindrical 
grid  with  a  constant  mesh  spacing  of  1  cm  in  both  the  radial  and  axial  directions.  Only  the  IMD  target 
plate  assembly  was  modeled  due  to  the  difficulty  of  modeling  the  guide  pistons  and  the  rest  of  the  IMD 
structure  in  CTH.  The  thickness  of  the  target  plate  assembly  was  adjusted  to  account  for  the  total  mass 
of  the  IMD  (5500  Ibm).  The  SESAME  EOS  for  4340  steel  was  used  for  the  target  plates.  Due  to  the  lack  of 
specific  values  for  4340  steel,  the  Johnson-Cook  default  values  for  steel  were  used  for  the  CTH 
constitutive  model.  The  JWL  EOS  was  used  to  model  the  explosive.  As  previously  stated,  the  soils  were 
modeled  using  the  HEP  model  fits  supplied  by  ERDC.  The  surrounding  air  was  modeled  using  the 
SESAME  tables  for  air.  Half  symmetry  was  used  to  reduce  the  problem  size.  Semi-infinite  boundary 
conditions  were  used  at  the  bottom,  top,  and  right  boundaries.  A  second  set  of  runs  was  repeated  using 
with  the  mesh  spacing  refined  to  5  mm  in  both  directions  to  quantify  the  effect  of  mesh  spacing  on  the 
results.  The  effects  of  CTH  inputs  such  as  the  multi-material  cell  pressure  formulation  (MMP)  and  cell 
fracture  (actually  spall)  parameters  were  also  investigated.  The  initial  material  configuration  for  all  the 
CTH  simulations  is  shown  Figure  2.13. 
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Figure  2. 13:  Initial  material  configuration  for  CTH  HEP  IMD  simulations 

Initial  simulations  were  run  using  the  default  CTH  MMP  option  (MMP0)  with  the  default  values  for 
material  strength  mixture  parameter  (MIX  1).  The  pressure  form  of  the  PFRAC  input  was  used  with  the 
PFRAC  value  for  the  soil  set  to  a  small  value  of  -l.E-4  dynes/cmA2.  In  the  absence  of  a  specific  material 
fracture  model  such  as  the  Johnson-Cook  model,  CTH  will  use  the  values  of  PFRAC  to  define  a  tensile 
pressure  at  which  a  material  is  considered  to  have  fractured  and  spalled.  CTH  uses  this  pressure  (or 
stress)  to  determine  when  void  material  is  introduced  into  the  cell.  This  induces  a  change  in  volume  in 
each  material  of  the  cell  that  has  a  pressure  more  negative  than  a  "relaxation"  value  that  is  a  fraction  of 
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the  fracture  criterion  for  the  cell.  CTH  also  allows  the  PFRAC  values  to  be  specified  as  a  principle  stress 
value  that  is  converted  to  a  pressure  value  inside  the  code.  However,  the  pressure  form  of  the  PFRAC 
input  is  the  preferred  input  format  for  CTH.  This  introduces  a  problem  with  the  HEP  model  due  to  the 
reference  of  the  pressures  generated  by  the  model  to  gage  pressure.  The  hydrostatic  (compression 
dependent)  term  in  the  HEP  EOS  allows  small  values  of  negative  pressure  (tension)  to  be  generated 
during  unloading.  However,  when  these  values  are  referenced  to  absolute  pressure  by  adding  the  gage 
pressure  (1  atm),  they  will  longer  be  negative.  The  results  obtained  to  date  indicate  that  the  simulations 
will  still  run  and  produce  realistic  results  but  the  exact  effect  of  using  the  pressure  based  PFRAC  is  still 
under  investigation.  Results  were  also  generated  using  the  stress  form  of  PFRAC  but  showed  only  small 
changes  in  results. 

2.4.4  Results  for  1  cm  mesh  with  strength 

The  impulse  of  the  IMD  target  obtained  for  the  three  soils  are  compared  with  the  experimental 
values  obtained  by  Ehrgot  [20]  and  the  EPIC  results  of  Moral  et  al.  [5]  for  the  default  6  mm  mesh  spacing 
given  is  shown  in  Figure  2.14.  On  examination  of  the  experimental  values  tabulated  in  the  Ehrgot's 
report,  it  was  found  that  the  experimental  value  for  peak  impulse  used  by  Moral  et  al.  for  the  dry  sand 
material  was  higher  (7740  N-s  versus  7273  N-s)  than  the  value  given  in  Ehrgot's  report.  The  lower 
experimental  value  for  impulse  of  the  dry  sand  is  shown  in  Figure  2.14.  For  the  dry  sand  material,  CTH 
HEP  under-predicts  the  impulse  by  about  15  per  cent.  EPIC  over-predicts  the  impulse  for  the  6  mm  mesh 
spacing.  For  the  wet  clay  material,  CTH  HEP  yields  a  better  correlation  with  experiment  than  EPIC. 
Impulse  for  the  wet  clay  material  is  slightly  over-predicted  by  CTH  HEP  and  under-predicted  by  EPIC.  For 
the  intermediate  soil  material,  both  EPIC  and  CTH  HEP  over-predict  the  impulse  in  comparison  with 
experiment  with  the  CTH  HEP  results  closer  to  experiment. . 
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Figure  2. 14:  Correlation  of  CTH  HEP  and  EPIC  impulse  on  IMD  for  course  mesh  with  experiment 
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Figure  2.15  is  taken  from  Reference  20  and  shows  the  measured  crater  dimensions  for  the  experiments. 
Material  plots  for  the  three  materials  at  15ms  simulation  time  are  given  in  Figures  2.16a-2.16c.  These 
results  demonstrate  the  importance  of  the  strength  and  failure  models  used  in  the  simulations  and  the 
effect  of  soil  cohesion  on  the  results.  The  wetter  materials  have  much  lower  shear  strength  but  higher 
cohesion  than  the  dry  sand.  This  leads  to  more  lateral  confinement  of  the  blast  for  a  longer  period  of 
time  for  the  wetter  materials  Visual  comparison  of  the  CTH  result  with  the  experimental  results  shown 
in  Figure  2.15  indicate  that  CTH  HEP  is  doing  a  good  job  of  predicting  the  crater  depth  and  as  well  as  the 
overall  shape  of  the  crater.  The  dry  sand  results  are  of  particular  interest  because  they  are  another 
indicator  of  the  effect  of  the  strength  model.  The  relatively  flat  section  of  the  crater  floor  is  captured 
along  with  the  shallow  slope  of  the  crater  wall.  However,  the  crater  diameter  appears  to  be  under¬ 
predicted.  Experience  has  shown  that  trying  to  interpolate  crater  geometry  from  CTH  material  plots 
should  be  done  with  caution  because  the  image  presented  can  be  greatly  affected  by  the  use  of  material 
discards  that  are  a  necessary  evil  in  many  CTH  calculations  to  insure  a  stable  solution.  A  better  method 
for  computing  the  crater  dimensions  based  on  material  volume  fraction  is  being  investigated.  The 
material  plots  for  the  wet  clay  and  intermediate  soil  materials  show  the  steeper  crater  slopes  and  lateral 
containment  of  the  blast  that  is  expected  due  to  their  higher  cohesions.  Crater  depths  for  these  two 
materials  are  also  close  to  the  experimental  values. 

2.4.5  Results  for  5  mm  mesh  with  strength 

In  order  to  determine  the  effect  of  mesh  density  on  results,  a  second  set  of  simulations  were  run 
with  the  CTH  mesh  size  decreased  to  5  mm  in  both  directions.  The  impulse  values  are  compared  with 
EPIC  results  using  a  3  mm  mesh  and  experiment  in  Figure  2.17.  The  effect  of  mesh  refinement  on  the 
CTH  results  is  to  shift  the  impulse  values  down  to  lower  peak  values.  This  moves  the  wet  clay  results 
slightly  below  the  experimental  values.  The  intermediate  soil  values  move  closer  to  the  experiment  but 
still  over-predict  the  impulse.  The  dry  sand  values  show  a  slight  downward  movement.  It  is  interesting  to 
note  that  the  EPIC  results  show  a  similar  effect  for  grid  refinement  from  6  mm  to  3  mm.  EPIC  continues 
to  provide  a  better  prediction  for  the  dry  sand  case  indicating  that  the  strength  and  failure  models  in 
EPIC  are  providing  a  more  realistic  result  for  this  soil.  In  contrast,  the  EPIC  and  CTH  results  for  the 
intermediate  soil  are  now  very  close. 


UNCLASSIFIED:  Distribution  Statement  A.  Approved  for  Public  Release 


30  |  P  a  g  e 


(wo)  A 


ENHANCEMENTS  AND  ANALYSIS  OF  CTH  SOFTWARE  FOR  UNDERBODY  BLAST 


Figure  2. 15:  Experimental  crater  dimensions  (Ref.  20) 
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Figure  2. 16:  Material  plots  at  15ms  for  three  IMD  soils  showing  crater  geometry 
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Figure  2. 17:  Correlation  of  CTH  HEP  and  EPIC  impulse  on  IMD  with  experiment  for  refined  mesh 


2.4.6  Effect  of  material  strength  on  impulse  and  crater  geometry:  Results  and  Discussion 

In  an  attempt  to  determine  the  cause  of  the  under-prediction  of  the  impulse  for  the  dry  sand  case,  a 
series  of  runs  were  made  in  which  the  HEP  strength  model  was  bypassed.  The  computed  values  of 
impulse  for  CTH  HEP  with  and  without  strength  are  shown  in  Figure  2.18.  There  is  a  slight  movement 
upward  toward  the  experimental  values  for  the  wet  soil.  The  intermediate  soil  results  move  upward  by  a 
large  amount.  For  the  dry  sand  case,  the  absence  of  material  strength  moves  the  CTH  results  closer  to 
the  EPIC  3  mm  predictions.  The  large  changes  in  the  intermediate  soil  and  dry  sand  results  indicate  the 
importance  of  the  material  strength  model  in  these  shallow  buried  simulations.  Another  example  of  this 
is  seen  in  the  material  plot  for  dry  sand  without  the  material  strength  shown  in  Figure  2.19.  The  crater 
dimensions  are  totally  erroneous.  The  dry  sand  results  confirm  that  a  better  implementation  of  the  HEP 
strength  model  is  needed  along  with  a  better  procedure  for  defining  material  failure  and  fracture.  It 
should  be  noted  that  these  results  were  generated  using  different  values  for  MMP  (MMP3)  and  PFRAC 
(principle  stress  form)  due  to  some  stability  problems  encountered  so  the  values  with  strength  will  differ 
slightly  from  the  previous  results.  The  effects  of  these  parameters  are  described  in  the  next  section. 
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Figure  2.  18:  Effect  of  material  strength  on  impulse  on  IMD  for  CTH  HEP  model  with  and  without  strength 
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Figure  2.  19:  Material  response  at  15  ms  of  IMD  DRY  SOIL  material  without  strength 


2.4.7  Effect  of  CTH  MMP  and  PFRAC  values:  Results  and  Discussion 

Figure  2.20  illustrates  the  effect  of  changing  the  CTH  MMP  and  PFRAC  parameters  on  the  computed 
impulse.  The  values  using  the  default  MMP  option  (MMP0)  with  a  pressure  based  PFRAC  value  of  -1.E4 
dynes/cmA2  are  compared  with  results  using  MMP3  and  a  stress  based  PFRAC  input  with  the  principle 
stress  for  fracture  assumed  to  be  -1.E4  dynes/cmA2.  The  CTH  MMP  options  determine  how  CTH 
distributes  volume  and  energy  for  multi-material  cells.  For  MMP0  this  distribution  is  a  linear  function  of 
volume  fraction.  MMP0  has  no  mechanism  for  relaxing  pressures  to  an  equilibrium  value.  By  contrast, 
MMP3  sets  a  limit  on  the  amount  of  PdV  work  that  is  allocated  based  on  a  cutoff  value  of  volume 
fraction.  In  addition,  pressure  relaxation  is  allowed  between  materials  in  a  cell.  It  can  be  seen  that  use  of 
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MMP3  and  a  stress  based  PFRAC  shifts  the  impulse  values  down  with  the  amount  of  shift  dependent  on 
the  material  strength.  This  result  indicates  that  more  testing  needs  to  be  done  to  quantify  the  best 
combination  of  MMP  and  PFRAC  for  CTH  HEP  runs. 


Effect  of  CTH  MMP  and  PFRAC  Models  on  Impulse 

Comparison  of  MMPO  and  Pressure  PFRAC  with  MMP3  and  Stress  PFRAC  Results 
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Figure  2.  20:  Correlation  of  impulse  on  IMD  for  CTH  HEP  model  with  different  values  of  PFRAC  and  MMP 


2.4.8  Effect  of  JC  Fracture  model  on  impulse  and  crater  dimensions:  Results  and  Discussion 

In  an  attempt  to  determine  the  cause  of  the  under-prediction  of  the  impulse  for  the  dry  sand  case,  a 
series  of  runs  were  made  to  evaluate  methods  to  improve  the  correlation  of  impulse  for  the  dry  sand 
material.  A  review  of  the  HEP  model  implementation  in  the  EPIC  code  revealed  that  EPIC  requires  an 
effective  strain  at  failure  as  an  input.  This  value  is  used  to  determine  when  EPIC  flags  a  material  as  failed 
and  then  bypasses  the  strength  calculation  for  that  material.  A  review  of  the  strength  and  failure  models 
in  CTH  showed  that  the  Johnson-Cook  fracture  model  could  provide  a  similar  capability  for  the  HEP 
model  if  only  the  first  of  the  five  standard  JC  fracture  model  inputs  (Dl)  was  set  to  a  small  strain  value 
and  the  spall  pressure  input  was  set  to  a  suitable  PRAC  value.  After  a  series  of  numerical  experiments,  it 
was  determined  that  a  Dl  value  of  4%  with  a  spall  pressure  (PF0)  of  -1.0E6  showed  substantial 
improvement  in  the  impulse  correlation  with  EPIC  and  the  IMD  experiments  at  the  cost  of  a  slightly 
worse  prediction  for  the  dry  sand  material.  For  the  other  two  IMD  materials,  using  the  JC  fracture  model 
had  either  little  or  no  effect  for  the  wet  clay  material  or  slightly  degraded  correlations  for  the 
intermediate  soil.  Figure  2.21  shows  the  impulse  correlation  using  the  JC  fracture  model.  Figure  2.22 
shows  the  predicted  crater  dimensions  for  the  three  soils.  These  results  illustrate  the  need  for  further 
research  into  an  appropriate  failure  algorithm  for  the  HEP  model. 
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Figure  2.  21:  Correlation  of  Impulse  on  IMD  for  CTH  HEP  and  EPIC  with  experiment  using  JC  fracture  model 
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Figure  2.  22:  Material  response  at  15ms  for  IMD  soils  using  CTH  HEP  with  JC  fracture  model 
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2.4.9  Simulation  of  the  ERDC  IMD  experiment  with  Bessette’s  IMD  model 

The  IMD  shallow  buried  explosion  test  for  the  WET  CLAY  material  was  rerun  using  Bessette's  more 
detailed  CTH  model  of  the  experiment.  Results  for  impulse  and  material  response  of  the  soil  generated 
by  the  CTH  HEP  model  were  compared  with  results  obtained  using  Bessette's  fit  for  the  IMD  WET  CLAY 
material  using  Kerley's  model.  The  CTH  material  plot  of  the  initial  material  configuration  of  the  detailed 
model  is  shown  in  Figure  2.23.  The  piston-plate  assembly  is  shown  in  gray.  The  darker  brown  material  is 
the  test  bed  soil  (WET  CLAY).  The  lighter  brown  material  is  the  in  situ  soil  that  is  modeled  as  a  silty  clay 
using  the  Kerley  model.  The  C4  charge  shown  in  red  was  modeled  using  the  History  Variable  Reactive 
Burn  model  (HVRB).  The  small  spot  of  green  material  is  a  PETN  booster  charge  that  is  used  to  initiate  the 
detonation  process  for  the  HVRB  model  and  is  modeled  using  the  JWL  EOS.  The  set  of  material  discards 
developed  by  Bessette  for  the  Kerley  model  was  used  for  both  tests.  The  surrounding  air  is  removed 
from  this  plot  for  clarity.  A  very  fine  2D  cylindrical  mesh  (2mm  cell  size)  is  used  in  the  areas  around  the 
explosive,  test  bed  soil,  and  IMD  piston-plate  assembly.  Mesh  stretching  is  used  to  extend  the 
boundaries  of  the  computational  domain  to  appropriate  distances  for  the  selected  boundary  conditions. 
The  total  number  of  computational  cells  used  in  the  simulation  totaled  over  7  million.  The  extent  of  the 
soil  materials  is  taken  to  be  semi-infinite.  All  runs  were  made  using  192  cores  on  the  MSU  Talon  cluster. 

Figure  2.24  compares  the  impulse  up  to  a  time  of  15  milliseconds  generated  by  both  the  CTH  HEP  and 
Bessette's  Kerley  model  fit  with  the  peak  impulse  reported  by  Ehrgott.  Both  models  are  seen  to  provide 
excellent  correlation  with  the  experiment.  Figures  25a-25f  show  the  material  response  for  both  models 
at  t=5.e-4,  2.e-3,  and  15e-3  seconds.  The  response  at  5e-4  and  2.e-3 
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Figure  2.  23:  Initial  material  configuration  for  Bessette's  model  of  IMD  experiment 


seconds  are  very  similar.  The  5e-4  second  shows  some  loss  of  material  at  the  top  of  the  soil  plug  being 
uplifted  by  the  explosion.  This  loss  is  most  likely  due  to  the  effects  of  the  material  discards  used  in  the 
problem.  The  late  time  material  plots  show  more  difference  than  the  early  time  plots  but,  as  stated 
previously;  care  must  be  taken  when  trying  to  interpret  information  on  crater  growth  from  CTH  material 
plots.  One  interesting  aspect  of  these  runs  was  the  dramatic  difference  in  run  times  required  for  each 
simulation.  As  seen  in  previous  tests,  the  CTH  HEP  implementation  was  dramatically  more  stable 
requiring  83,311  cycles  and  30,320  cpu  seconds  to  reach  completion  compared  to  229,492  cycles  and 
82,090  cpu  seconds  for  the  Kerley  model.  This  difference  in  stability  is  further  illustrated  in  Figure  2.26 
that  compares  the  stable  time  step  computed  by  CTH  for  both  models. 

2.4.10  Summary  of  ERDC  IMD  simulations:  Results  and  Discussion 

The  results  from  the  ERDC  IMD  simulations  revealed  that  CTH  HEP  provides  excellent  to  good 
correlations  with  experiment  and  EPIC  calculations  for  the  wetter  materials  with  low  shear  strength.  The 
correlations  are  not  as  good  as  EPIC  for  the  higher  strength  dry  sand  material.  These  results  indicate  that 
more  work  needs  to  be  done  on  the  CTH  HEP  strength  model  implementation  and  material  failure 
models.  Several  options  are  being  explored  including  using  a  CTH  "modern"  strength  model  that 
computes  the  deviatoric  stresses  for  individual  materials  and  implementing  an  effective  strain  failure 
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criterion  similar  to  the  one  used  in  EPIC.  However,  it  is  felt  that  the  correlations  obtained  for  all  three 
materials  are  sufficiently  accurate  for  use  in  realistic  simulations  of  shallow  buried  explosions. 

IMD  Impulse  Versus  Time 


Correlation  of  Computed  Impulse  with  Experiment  for  Bessette  s  IMD  Model 


Figure  2.  24:  Comparison  of  Computed  IMD  Impulse  with  Experiment  for  HEP  and  Bessette's  fit 
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2.25a)  Bessette's  fit  (t=5e-4s) 
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IMD  Wet  Clay  Test  -  T  =  2.00e-03  Seconds 

400  I  ■  ■  ■  -  i  . 

300  1 
200  1 
100  1 


-300  -200  -100  0  100  200  300 

2.25c)  Bessette's  fit  (t=2e-3s) 
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2.25d)  HEP  (t=2e-3s) 


IMD  Wet  Clay  Test  -  T  =  1 .50e-02  Seconds 
400  . 

300  1 
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2.25e)  Bessette's  fit  (t=15e-3s) 

Figure  2.  25:  Comparison  of  material  plots  for  Bessette's  fit  and  HEP  at  three  times 
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CTH  Stable  Time  Step  Versus  Simulation  Time 

Comparison  of  CTH  Stable  Time  Step  for  HEP  and  Bessette  Sesame  IMD  WET  CLAY  Fit 
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Figure  2.  26:  Comparison  of  stable  time  step  for  HEP  and  Bessette  WET  CLAY  runs 
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Section  3:  Development  of  a  one-way  CTH  to  LS-DYNA  coupling 
procedure 

Over  the  years  several  researchers  have  attempted  to  provide  a  means  of  extracting  results  from  CTH 
calculations  for  use  in  Lagrangian  Finite  Element  structural  dynamics  codes  such  as  LS-DYNA  and 
DYNA3D.  Of  these  the  most  successful  but  also  the  most  complex,  is  the  ZAPOTEC  code  from  Sandia  [32] 
that  links  CTH  to  the  Sandia  PR0NT03D  structural  dynamics  code.  In  ZAPOTEC,  the  coupling  is  referred 
to  as  two-way,  meaning  that  both  codes  run  as  a  unit  with  information  such  as  material  displacement, 
movement  of  mesh  nodes,  and  forces  and  moments  are  exchanged  at  either  every  time  step  or  some 
multiple  of  the  stable  time  step  required  by  each  code  component.  The  two-way  approach  provides  the 
most  physically  realistic  simulations  of  Fluid  Structure  Interaction  at  a  cost  of  code  complexity  and  run 
time.  An  alternate  approach  is  commonly  referred  to  as  a  one-way  coupling  in  which  each  code  is  run  in 
separate  steps  with  information  from  one  code  saved  in  a  format  the  other  code  can  read.  The 
underlying  assumption  of  the  one-way  coupling  is  that  the  response  of  the  structure  is  much  slower 
than  the  time  scale  of  the  blast  or  other  dynamic  loading  event  modeled  by  CTH.  Therefore,  it  is 
assumed  that  extracting  pressures  from  codes  such  as  CTH  up  to  a  time  shorter  than  the  simulation  time 
required  to  model  the  complete  structural  response  and  then  using  those  pressures  as  time  dependent 
loads  on  the  structure  will  provide  an  acceptable  approximation  to  the  total  force  and  impulse  imparted 
to  the  structure  at  a  reduce  cost  when  compared  with  the  two-way  approach.  Another  underlying 
assumption  that  is  common  to  most  one-way  couplers  is  that  the  material  to  be  loaded  can  be 
considered  a  rigid  material  in  the  code  (CTH  in  this  case)  that  generates  the  pressures  used  by  the 
structural  dynamics  code.  The  one-way  coupling  approach  is  the  one  selected  for  the  current  work 
although  it  might  be  more  accurate  to  refer  to  the  procedure  developed  as  a  load-curve  extractor.  The 
following  sections  describe  the  implementation  of  a  one-way  coupler  between  CTH  and  the  LS-DYNA 
structural  dynamics  code  and  the  testing  done  to  verify  the  method. 

3.1  Approach 

The  approach  taken  in  the  current  CTH-LS-DYNA  coupler  takes  its  inspiration  (and  some  of  its  code) 
from  a  previous  one-way  coupler  developed  by  Dr.  David  Littlefield  (UAB)  [27],  Littlefield's  approach 
allowed  users  to  define  the  information  required  to  generate  load  curves  for  LS-DYNA  in  standard 
KEYWORD  format.  Software  embedded  in  CTH  would  then  read  the  LS-DYNA  input  and  generate  an 
output  file  containing  for  load  curves  (DEFINE_CURVE  input)  in  a  format  that  could  be  pasted  directly 
into  the  LS-DYNA  inputs  without  modification. 

Initial  attempts  to  implement  the  current  coupler  tried  to  use  the  software  developed  by  Littlefield 
for  an  earlier  version  of  CTH  but  these  proved  unsuccessful  due  to  differences  in  CTH  9.1  and  the 
previous  version.  The  primary  difference  was  that  the  rigid  material  model  developed  by  Littlefield  that 
was  a  basic  component  of  his  coupler  could  not  be  made  to  work  in  CTH  9.1.  Therefore,  it  was  decided 
to  keep  some  of  the  basic  components  of  Littlefield's  coupler  but  replace  his  rigid  material  model  with 
the  native  rigid  CTH  rigid  material  model.  An  initial  implementation  was  developed  that  assumed  the 
elements  in  the  LS-DYNA  model  were  solid  elements  since  Littlefield's  coupler  only  supported  solid 
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elements.  As  work  on  this  task  progressed,  it  became  apparent  that  Littlefield's  software  would  have  to 
be  extended  to  allow  users  to  use  shell  elements  in  their  LS-DYNA  models.  This  requirement  led  to  a 
complete  rewrite  of  most  of  the  existing  Littlefield  software.  A  description  of  the  new  and  modified 
routines  that  make  up  the  coupler  is  given  in  Appendix  B. 

A  second  major  difference  in  the  current  implementation  and  the  Littlefield  implementation  is  the 
use  of  a  time-average  pressure  based  on  an  approximation  of  the  impulse  delivered  to  the  rigid  target 
material  instead  of  the  instantaneous  pressure.  The  initial  implementation  of  the  coupler  used  the 
instantaneous  pressure  approach  and  it  appears  to  yield  satisfactory  results  for  many  cases.  However, 
the  time  average  approach  is  the  one  used  in  a  one-way  coupler  between  CTH  and  the  PRESTO 
structural  dynamics  component  of  the  Sandia  SIERRA  framework  [28],  A  careful  review  of  the  inner 
workings  of  the  CTH  rigid  material  model  revealed  that  the  time-averaged  pressure  approach  is  the  one 
used  in  CTH  to  apply  forces  on  a  rigid  material.  Subsequent  testing  of  the  time-averaged  approach  using 
the  DRDC  plate  experiment  to  be  described  later  showed  improved  results  over  the  instantaneous 
pressure  approach. 

The  basic  steps  required  for  a  coupled  analysis  are  defined  in  the  next  section. 

3.1.1  Steps  required  for  a  typical  CTH-LS-DYNA  coupled  analysis 

A  coupled  analysis  using  CTH  and  LS-DYNA  typically  involves  the  following  steps: 

1.  An  LS-DYNA  model  is  generated  and  the  surfaces  to  be  loaded  by  CTH  are  defined  by  LS-DYNA 
LOAD_SEGMENT  or  LOAD_SHELL  inputs.  An  LS-DYNA  keyword  input  file  containing  the  LOAD_SEGMENT 
or  LOAD_SHELL  input  along  with  the  NODE  and  ELEMENT  input  that  correspond  to  the  LOAD  input 
needs  to  be  created  to  use  as  input  to  the  CTH  coupler 

2.  If  the  geometry  of  the  LS-DYNA  model  cannot  be  replicated  using  the  native  CTH  material  insertion 
facility  (called  diatom),  the  LS-DYNA  model  should  be  exported  using  a  preprocessor  such  as  the  LS- 
DYNA  LS-PREPOST  as  a  Stereo  Lithography  (STL)  format  file.  Note  the  geometry  exported  only  needs  to 
be  the  surfaces  that  will  be  exposed  to  the  CTH  generated  loads. 

3.  A  CTH  input  file  is  created  for  the  loading  side  of  the  problem.  It  should  contain  a  representation  of 
the  target  geometry  and  inputs  in  a  new  CTH  input  block  (delimited  by  a  fern  and  endfem  statement)  to 
define  the  input  LS-DYNA  keyword  file,  output  file  that  will  contain  the  appropriate  DEFINE_CURVE 
corresponding  to  each  LOAD_SEGMENT  or  LOAD_SHELL  input,  and  other  inputs  such  as  the  frequency 
that  the  load  curves  are  dumped  to  the  output  file. 

4.  CTH  is  run  for  a  predefined  simulation  time  to  generate  the  load  curve  output  file.  This  load  curve  file 
is  then  pasted  into  the  LS-DYNA  input  files  (usually  as  an  LS-DYNA  *INCLUDE  file) 

5.  The  final  step  is  to  run  the  LS-DYNA  simulation 

A  more  detailed  description  of  the  coupler  is  given  in  the  User's  Manual  that  accompanies  the  release 
of  the  enhanced  version  of  CTH  9.1. 
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3.1.2  Initial  implementation  tests 

Some  initial  test  of  the  CTH  component  of  the  coupler  were  run  to  verify  the  software  that  read  the 
LS-DYNA  input  files,  to  examine  the  accuracy  of  the  procedure  used  to  extract  the  instantaneous 
pressures  used  in  the  initial  coupler  implementation  and  to  test  the  parallel  implementation  of  the 
software.  The  implementation  of  the  new  coupling  software  in  CTH  9.1  was  initially  verified  using  two 
test  cases:  an  explosively  loaded  rigid  witness  plate  and  an  example  LSDYNA  case  that  models  blast 
loads  on  a  cartoon  "car".  The  rigid  plate  case  was  run  in  both  serial  and  parallel  modes.  The  load-blast- 
car  case  was  run  in  parallel  only  due  to  the  size  of  the  grid  required  to  contain  the  car  geometry.  Both 
cases  are  described  in  the  following  sections.  The  results  presented  used  the  instantaneous  pressure 
instead  of  the  time-  averaged  pressures  used  in  the  final  implementation  of  the  coupler. 

3. 1.2.2  Rigid  plate  case:  Results  and  Discussion 

The  modifications  to  CTH  9.1  were  first  verified  using  a  simple  test  case  that  consisted  of  an 
explosively  loaded  plate.  The  plate  is  10  cm  high,  10  cm  wide  and  1  cm  thick.  The  plate  is  loaded  by  a 
spherical  charge  of  TNT  with  a  radius  of  7  cm.  This  results  in  a  charge  mass  of  about  2.3  kg  which  is 
equivalent  to  a  charge  weight  of  about  5  Ibf.  The  standoff  distance  of  the  plate  from  the  center  of  the 
TNT  sphere  is  10  cm.  The  CTH  programmed  burn  model  and  the  JWL  equation  of  state  (EOS)  for  TNT  are 
used  to  model  the  explosive.  The  surrounding  air  is  modeled  using  the  SESAME  tabular  EOS  for  air.  The 
plate  is  considered  a  rigid  material.  The  CTH  mesh  contained  80x80x80  cells  with  a  constant  0.25  cm 
mesh  spacing  in  all  directions.  The  Z=0  and  X=0  planes  are  specified  to  be  a  reflective  boundaries.  The 
other  boundaries  are  specified  to  represent  a  semi-infinite  medium.  The  simulation  time  for  this  case 
was  3.E-4  seconds.  This  test  case  was  run  in  both  serial  and  parallel  mode  on  a  Linux  workstation  with  6 
core  AMD  X86-64  processor.  CTH  was  compiled  with  version  11  of  the  Intel  Fortran  and  C  compilers. 
OpenMPI  was  used  for  the  MPI  libraries.  The  parallel  runs  used  6  cores.  A  material  plot  showing  the 
initial  configuration  of  the  materials  is  shown  in  Figure  3.1.  A  contour  plot  of  pressure  at  2.e-5  seconds  is 
shown  in  Figure  3.2.  The  spherical  blast  front  is  seen  to  impact  along  the  length  of  the  plate  with  the 
outer  parts  of  the  front  wrapping  around  the  sides  of  the  plate  and  reflecting  on  the  symmetry  plane. 
The  plate  is  clearly  seen  to  be  non-responding  indicating  that  CTH  is  treating  the  plate  material  as  rigid. 

The  LS-DYNA  LSPREPOST  tool  was  used  to  create  an  initial  LS-DYNA  model  for  the  plate.  The  LS-DYNA 
mesh  for  the  plate  consists  of  10x10x10  hexahedral  solid  elements  with  1331  nodes.  Each  cell  face  on 
the  Z=10  surface  is  specified  to  be  a  LOAD_SEGMENT.  The  keyword  input  file  generated  by  LSPREPOST 
for  the  plate  model  was  used  as  the  input  file  for  the  coupling  algorithm.  Therefore,  there  will  be  one 
hundred  load  curves  in  the  resulting  coupling  output  file.  The  pressure  time  histories  extracted  by  the 
coupling  procedure  were  compared  with  CTH  tracer  data  for  the  tracer  locations  that  corresponded  to 
the  LOAD_SEGMENT  surface  centroids  as  well  as  tracers  at  cell  centers  in  the  cell  closest  to  the 
LOAD_SEGMENT  location  at  the  LOAD_SEGMENT  locations  shown  in  Figure  3.3.  For  the  given  CTH  mesh 
configuration,  the  LOAD_SEGMENT  surface  (Z=10)  lies  on  a  cell  face.  By  default,  CTH  will  try  to 
interpolate  the  cell  center  data  to  the  cell  face.  For  the  rigid  material,  this  lead  to  the  discrepancy  in  the 
extracted  load  curve  and  the  interpolated  data  shown  in  Figure  3.4.  The  extracted  curve  peak  is  almost 
exactly  twice  the  peak  for  interpolated  tracer  data.  In  trying  to  reconcile  this  discrepancy,  it  was 
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determine  that  CTH  sets  the  pressure  inside  a  rigid  material  to  zero.  Therefore  ,  a  simple  average  of  the 
cell  values  on  either  side  of  the  rigid  surface  would  lead  to  an  interpolated  value  of  exactly  one  half  of 
the  value  of  pressure  at  the  first  cell  center  below  the  surface.  The  almost  exact  correlation  between 
the  extracted  data  and  tracer  data  at  the  corresponding  cell  center  shows  that  the  pressure  extraction 
method  for  the  rigid  material  is  using  the  first  cell  center  off  of  the  surface.  Correlations  at  the  other 
three  LOAD_SEGMENT  locations  shown  in  Figure  3.5a-3.5c  verify  this.  Figure  3.6  compares  the 
extracted  load  curve  data  for  the  serial  and  parallel  runs  at  the  center  of  the  plate.  The  correlation 
indicates  that  the  parallel  version  of  the  code  is  generating  the  same  results  as  the  serial  version. 


Materials  Time  0.00e+00  Seconds 


Figure  3. 1:  Initial  material  configuration  for  rigid  plate  test 
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Pressure  at  2.00e-05  Seconds 


Pressure 


Figure  3.  2:  Rigid  plate  pressure  contours  at  t=2.e-5  seconds 


■  Geom  Parts 


LS-DYNA  keyword  deck  by  LS-PrePost 


Figure  3.  3:  LOAD_SEGMENT  locations  used  in  pressure  correlations 
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Rigid  Plate  Pressure  vs  Time 


Load  Segment  Location  :  X=0.5  cm  Y=0.5  cm  Z=10  cm 


Figure  3.  4:  Comparison  of  interpolated  and  extracted  pressures  with  cell  center  data 


Rigid  Plate  Pressure  vs  Time 


Load  Segment  Location  :  X=0.5  cm,  Y=4.5  cm,  Z=10  cm 


3.5a)  Pressures  at  X=0.5,  Y=4.5 
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Rigid  Plate  Pressure  vs  Time 

Load  Segment  Location  :  X=4.5  cin,  Y=4.5  cm,  Z=10  cm 


3.5b)  Pressures  at  X=4.5,  Y=4.5 


Rigid  Plate  Pressure  vs  Time 


Load  Segment  Location  :  X=9.5  cm,  Y=4.5  cm,  Z=10  cm 


3.5a)  Pressures  at  X=9.5,  Y=4.5 

Figure  3.  5:  Comparison  of  extracted  pressures  will  cell  center  pressures  at  3  X  locations 
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Comparison  of  Serial  and  Parallel  Extracted  Load  Curves 


Load  Segment  Location  :  X=4.5  cin,  Y=4.5  cm,  Z=10  cm 


Figure  3.  6:  Comparison  of  serial  and  parallel  simulation  pressures 


3. 1.2. 3  LOAD-BLAST-CAR  case:  Results  and  Discussion 

In  order  to  validate  the  coupling  procedure  for  a  large  scale  problem  and  to  verify  that  STL  input  can 
be  used  to  define  a  rigid  material  in  CTH  9.1,  a  test  case  based  on  the  "cartoon"  automobile  geometry 
described  in  the  first  quarter  report  was  constructed.  The  geometry  and  corresponding  LS-DYNA 
keyword  input  file  was  taken  from  a  standard  LS-DYNA  example  problem  (load-blast-car.k)  used  to 
illustrate  the  use  of  the  LS-DYNA  LOAD_BLAST  facility.  The  input  was  modified  to  express  the  geometry 
in  CTH  units  (cm  versus  m)  and  to  swap  the  Y  and  Z  axis  to  be  consistent  with  the  preferred  coordinate 
system  for  3D  calculations  in  CTH  (Z  axis  is  up  and  Y  is  spanwise).  The  LSDYNA  generated  STL  file  is 
shown  in  Figure  3.7  for  the  initial  geometry  with  the  "tires"  removed.  The  modified  geometry  was  used 
in  LSPREPOST  to  generate  LOAD_SEGMENT  input  at  the  36  face  segments  shown  in  Figure  3.8a. 

A  CTH  input  deck  was  constructed  that  consisted  of  a  large  mesh  (900cm  x  400  cm  x  524  cm)  with  2 
cm  mesh  spacing  in  all  directions.  This  leads  to  a  mesh  with  450  x  200  x  262  interior  cells.  The  Z=-24 
plane  is  taken  as  a  reflecting  boundary.  All  other  exterior  boundaries  are  set  to  model  a  semi-infinite 
domain.  Although  the  mesh  is  large,  the  2  cm  spacing  is  not  fine  enough  in  the  region  of  the  explosive  to 
capture  the  generated  blast  wave  accurately.  However,  it  is  sufficient  for  the  objectives  of  this  case; 
mainly,  to  verify  the  coupling  algorithm  for  a  large  scale  case  and  demonstrate  the  use  of  the  STL 
geometry  from  LSPREPOST  to  define  the  rigid  material  in  CTH.  The  explosive  charge  consisted  of  a 
sphere  of  about  100  kg  of  TNT.  The  diameter  of  the  sphere  is  48  cm.  The  detonation  point  of  the 
explosive  was  at  X=340,  Y=100,  Z=-24.0  which  is  just  in  front  of  the  left  front  "tire"  of  the  car.  The 
explosive  was  modeled  using  the  JWL  equation  of  state.  Simulations  were  run  on  in  parallel  using  60 
processors  on  the  MSU  HPCC  Talon  cluster.  Pressure  time  histories  were  saved  and  extracted  every 
microsecond. 
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A  material  plot  for  a  2D  sections  of  the  grid  in  the  Y=0  cm  symmetry  plane  at  2.8e-4  seconds  is  shown 
in  Figures  3.9.  This  plot  illustrates  the  ability  of  the  CTH  STL  input  option  to  preserve  the  key  features  of 
the  geometry  when  it  is  inserted  into  a  3D  Cartesian  mesh. 

Results  comparing  extracted  pressure  time  histories  versus  data  at  CTH  tracer  locations  are 
presented  for  the  six  stations  shown  in  Figure  3.8b.  These  results  are  shown  in  Figure  3.10a-3.10e. 
Figure  3.10a  again  shows  the  discrepancy  between  interpolated  and  cell  center  data.  However,  the 
other  plots  in  Figure  3.9  show  a  near  exact  correlation  between  the  extracted  values  and  the  values  for 
tracers  defined  at  the  cell  center  closest  to  the  corresponding  LOAD_SEGMENT  location. 

These  results  indicate  that  the  coupling  algorithm  and  STL  geometry  definition  procedure  works  for  a 
realistic  large-scale  simulation.  In  particular,  the  baseline  coupling  algorithm  is  shown  to  work  with  the 
default  CTH  rigid  material  model. 


Figure  3.  7:  STL  Representation  of  cartoon  car 
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3.8a)  LOAD_SEGMENT  elements  defined  for  cartoon  car 


3.8b)  LOAD_SEGMENT  locations  used  in  pressure  correlations 

Figure  3.  8:  LOAD_SEGMENT  elements 
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Materials  Time  2.80e-04  Seconds 


Figure  3.  9:  Material  plot  in  Y=0  plane  at  2.8e-4  seconds 


Extracted  and  Interpolated  CTH  Tracer  Pressure  vs  Time 

Load  Segment  Location  :  X=333.  cm.  Y=95.  cm,  Z=0.  cm 
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3.10a)  Pressures  at  X=333,  Y=95,  Z=0 
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Extracted  and  CTH  Tracer  Pressure  vs  Time 


Load  Segment  Location  :  X=320  cm,  Y=-53.33  cm,  Z=30.  cm 


3.10b)  Pressures  at  X=320,  Y=-53,  Z=30 


Extracted  and  CTH  Tracer  Pressure  vs  Time 


Time  (sec) 


3.10c)  Pressures  at  X=357.5,  Y=0,  Z=30 
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Extracted  and  CTH  Tracer  Pressure  vs  Time 

Load  Segment  Location  :  X =357.6  cm,  Y=84.6  cm,  Z=30. 


3.10d)  Pressures  at  X=357.6,  Y=84.6,  Z=30 


Extracted  and  CTH  Tracer  Pressure  vs  Time 

Load  Segment  Location  :  X =427.5  cm,Y=-53.33  cm,  Z=30  cm 


3.10e)  Pressures  at  X=427.5,  Y=-53.33,  Z=30 
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Extracted  and  CTH  Tracer  Pressure  vs  Time 
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3.10f)  Pressures  at  X=457.54,  Y=89.1,  Z=30 
Figure  3. 10:  Comparison  of  extracted  pressures  and  cell  center  pressures  at  6  locations 


3.2  Full  scale  tests  of  the  one-way  coupler 

A  full  scale  test  of  the  coupler  was  delayed  by  implementation  of  the  HEP  model.  When  testing 
resumed,  it  focused  on  two  test  cases.  At  the  suggestion  of  TARDEC/CASSI  Analytics,  the  first  full  scale 
test  of  the  CTH  -  LS-DYNA  coupler  attempted  to  simulate  results  for  the  DRDC  plate  experiment  [24].  A 
second  full  scale  test  attempted  to  model  the  Generic  Hull  experiment  used  in  validating  the  Loci- 
Blast/LS-DYNA  coupler  developed  in  SimBRS  WD0034  [29].  The  DRDC  plate  experiment  has  been 
subjected  to  a  variety  of  analyses  over  the  last  several  years  including  use  of  LS-DYNA  ALE,  LS-DYNA  SPH 
and  FEM  coupling,  and  PAM-SHOCK  SPH-FEM  coupling  [30,  31].  The  PAM-SHOCK  SPH-FEM  approach 
was  also  applied  to  the  Generic  Hull  experiment  in  Reference  30.  The  next  sections  describe  results  from 
the  using  CTH  alone  and  the  CTH-LS-DYNA  coupler  to  simulate  the  DRDC  plate  experiment  and  the  CTH- 
LS-DYNA  coupler  to  simulate  the  Generic  Hull  experiment. 

3.2.1  Simulation  of  the  DRDC  plate  experiment:  Results  and  Discussion 

LS-DYNA3D  ALE  and  LOAD_BLAST  models  of  the  target  plate  and  support  structure  geometries  were 
provided  by  TARDEC  and  used  to  develop  a  CTH  input  deck.  On  reviewing  the  LS-DYNA  input  decks,  it 
was  found  that  shell  elements  were  used  in  the  LS-DYNA  analysis.  This  complicated  the  analysis  since 
CTH  9.1  does  not  handle  shell  elements.  Plate  thicknesses  were  added  to  the  LS-DYNA  shell  element 
definitions  to  generate  a  3D  model  of  the  materials  inserted  into  the  CTH  mesh.  For  these  simulations, 
the  top  plate  was  ignored  along  with  the  side  plates  used  to  constrain  the  lateral  motion  of  the  support 
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frame.  The  soil  was  modeled  using  the  CTH  HEP  model  for  dry  sand.  The  dry  sand  material  density  does 
not  match  the  density  of  the  soil  used  in  the  LS-DYNA  simulations  but  was  chosen  based  on  the  results 
from  the  ERDC  IMD  experiment  simulations  and  represents  a  lower  bound  on  impulse  that  will  be 
applied  by  the  soil  and  explosion  to  the  target  plate.  The  target  plate  used  was  the  5083  aluminum  plate 
used  in  the  LS-DYNA  inputs  we  received  from  TARDEC.  The  plate  thickness  used  was  the  experimental 
value  of  3.175  cm.  As  reported  in  Reference  24,  the  explosive  charge  consisted  of  a  6  kg  cylinder  of  C4 
with  a  diameter  of  25.4  cm  and  a  thickness  of  7.62  cm.  The  DOB  is  taken  to  be  to  the  top  of  the  charge 
and  was  set  at  5  cm.  The  default  Johnson-Cook  constitutive  model  for  5083  aluminum  was  used  with 
the  individual  parameters  adjusted  to  match  the  LS-DYNA  input  values.  The  Johnson-Cook  Fracture 
model  was  also  used.  The  SESEME  Equation  of  State  values  for  2024  aluminum  with  density  adjusted  to 
match  the  5083  value  was  used  for  the  target  plate  due  to  the  lack  of  a  default  EOS  model  for  5083 
aluminum.  The  SESAME  EOS  tables  for  4340  steel  were  used  for  the  EOS  input  and  the  default  Johnson- 
Cook  steel  values  were  used  for  the  constitutive  model.  The  surrounding  air  was  modeled  using  the 
SESAME  EOS.  Quarter  symmetry  was  used  to  reduce  the  problem  size.  However,  the  thin  plates  used  in 
the  geometry  require  a  fine  mesh  to  resolve  the  material.  Therefore  a  graded  mesh  using  5  mm  mesh 
spacing  in  the  immediate  region  of  the  test  geometry  and  the  charge  was  used  with  mesh  stretching 
used  to  extend  the  mesh  to  appropriate  outer  boundaries.  Even  with  a  mesh  grading,  the  resulting  mesh 
contains  280x280x480  (37,632,000)  cells.  Due  to  the  size  of  the  grid,  these  simulations  were  run  on  the 
MSU  HPCC  Talon  cluster  using  192  cores.  A  diagram  of  the  DRDC  plate  experiment  test  configuration 
taken  from  Reference  24  is  shown  in  Figure  3.11.  A  3D  view  of  the  initial  CTH  configuration  is  shown  in 
Figure  3.12.  An  additional  1  metric  ton  of  weight  that  sits  on  top  of  the  frame  is  not  shown.  An  initial  set 
of  CTH  alone  simulations  were  run  prior  to  running  the  CTH-LS-DYNA  coupler  to  examine  the  effect  of 
mesh  spacing  on  run  times,  the  stability  of  the  simulation  using  the  HEP  soil  model,  and  to  see  how 
accurately  CTH  could  model  the  DRDC  plate  response  on  its  own. 


Figure  3. 11:  DRDC  plate  experiment  test  configuration  (Reference  24) 


3.2.2  LS-DYNA  ALE  and  CTH  alone  simulation  of  DRDC  plate  experiment 

Prior  to  running  the  CTH  simulations,  the  LS-DYNA  input  decks  obtained  from  TARDEC  were  run  using 
both  the  ALE  and  CONWEP  inputs.  Figure  3.13  shows  the  plate  deformation  obtained  using  the  CONWEP 
based  LOAD_BLAST  facility  in  LS-DYNA  to  provide  the  blast  loads  to  the  structure.  For  the  CONWEP  case 
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the  plate  is  seen  to  deform  but  not  fail  which  matches  the  experiment.  Next  the  LS-DYNA  ALE  case  was 
run  using  the  default  material  properties  for  the  target  plate  obtained  from  TARDEC.  With  these  inputs, 
the  plate  is  seen  in  Figure  3.14  to  exhibit  fracture  under  the  blast  load  that  is  not  indicated  in  the 
experiment.  Figure  3.15  shows  the  response  obtained  from  running  CTH  alone  using  the  same  Johnson- 
Cook  inputs  for  the  plate  used  in  the  LS-DYNA  ALE  case.  The  CTH  results  show  an  even  more  dramatic 
plate  fracture  response  than  LS-DYNA.  Several  numerical  experiments  were  conducted  to  determine  the 
cause  of  the  fracture.  These  experiments  indicated  that  the  parameters  used  in  the  Johnson-Cook 
Fracture  model  were  resulting  in  premature  failure  of  the  plate.  These  parameters  were  modified  to  set 
the  leading  value  (Dl)  to  150%  and  the  rest  of  the  parameters  to  zero.  This  modification  was  similar  to 
that  used  with  the  JC  fracture  model  application  to  the  HEP  soil  model.  The  effect  of  this  modification  is 
to  turn  off  all  rate  and  thermal  effects  in  the  fracture  model.  Although  there  is  no  physical  basis  for  this 
modification,  the  results  presented  in  Figures  3.16-3.18  show  that  the  change  in  fracture  model  inputs 
leads  to  improved  prediction  of  the  deformation  response  for  both  LS-DYNA  and  CTH.  The  material 
deformation  responses  show  in  Figures  3.16  and  3.17  are  now  consistent  with  the  LS-DYNA  CONWEP 
results.  The  plate  centerline  displacement  time  history  is  show  in  Figure  3.18  for  both  the  LS-DYNA  ALE 
and  CTH  HEP  simulations.  The  predicted  value  of  displacement  at  6  ms  compares  favorably  with  the  300 
mm  value  measured  in  the  experiment.  These  results  illustrate  the  importance  of  the  fracture  model  in 
predicting  plate  response.  Although  the  IMD  DRY  SAND  material  yielded  surprisingly  good  results  for  the 
centerline  displacement  for  the  CTH  alone  tests,  this  should  not  be  considered  an  indication  that  it  is  an 
appropriate  analog  for  the  soil  used  in  the  DRDC  experiments. 


Time  O.OOe+OO  Seconds 


Figure  3.  12:  Initial  CTH  material  configuration  for  DRDC  plate  model 
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Figure  3. 13:  DRDC  Plate  deformation  for  LS-DYNA  LOAD_BLAST  analysis 


Figure  3. 14:  DRDC  Plate  deformation  for  LS-DYNA  ALE  analysis  with  default  JC  fracture  values 
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Figure  3. 15:  DRDC  plate  deformation  for  CTH  HEP  analysis  with  default  JC  fracture  values 
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Figure  3. 16:  DRDC  plate  deformation  for  LS-DYNA  ALE  analysis  with  modified  JC  fracture  values 
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Figure  3. 17:  DRDC  plate  deformation  for  CTH  HEP  analysis  with  modified  JC  fracture  values 


DRDC  Plate  Centerline  Displacment  vs  Time 


Comparison  of  CTH  HEP  and  LS-DYNA  ALE  Results  with  Modified  Plate  JC  Fracture  Inputs 


Figure  3.  18:  Comparison  of  predicted  DRDC  plate  centerline  displacement  versus  time 
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3.2.3  CTH  to  LS-DYNA  one-way  coupling  results  for  DRDC  plate  experiment 

Two  sets  of  simulations  of  the  DRDC  plate  experiment  were  performed  using  the  CTH/LS-DYNA 
coupler.  In  the  first  set  of  simulations,  the  instantaneous  pressures  on  the  plate  were  used  to  define  the 
load  curves.  In  the  second  set  of  simulations,  the  time  average  approach  to  defining  the  load  curve 
pressures  was  used.  The  CTH  inputs  for  the  DRDC  geometry  developed  for  the  CTH  alone  calculations 
were  used  with  the  exception  that  the  supporting  legs  were  removed  to  match  the  definition  of  the  load 
curves  generated  by  LS-PREPOST.  The  5  mm  mesh  spacing  used  in  the  CTH  alone  case  was  repeated  in 
the  coupled  cases.  CTH  simulations  were  first  run  to  5e-3  seconds  and  the  extracted  load  curves  were 
fed  into  an  LS-DYNA  simulation  that  ran  for  8e-3  seconds.  Figure  3.19  compares  the  displacement  of  the 
plate  at  element  locations  in  the  symmetry  plane  of  the  plate  from  the  plate  center  to  the  plate  edge 
computed  using  both  the  instantaneous  and  time  averaged  CTH  pressures  with  experiment.  It  is  seen 
that  the  time  average  results  gives  a  much  better  correlation  with  experiment  for  this  case.  The 
overprediction  of  the  deformation  near  the  plate  edge  with  the  time-averaged  pressures  needs  to  be 
investigated  further.  A  contour  plot  of  displacement  is  shown  in  Figure  3.20. 


DRDC  Plate  Vertical  Displacement  Vs  Distance  from  Plate  Center 

Correlation  of  CTH  Pressure  Load  Curve  Response  with  Experiment  (0.5cm  CTH  mesh) 


Distance  from  Plate  Center  (mm) 


Figure  3. 19:  Comparison  of  plate  deflection  using  time-averaged  and  instantaneous  CTH  pressures  with  experiment 
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Figure  3.  20:  Displacement  contours  using  time-averaged  CTH  pressures 


3.2.4  Simulation  of  Generic  Hull  experiment  using  CTH  to  LS-DYNA  one-way  coupling: 

Results  and  Discussion 

The  TARDEC  Generic  Hull  experiment  was  designed  to  provide  a  notional  geometry  for  underbody 
blast  analysis  and  to  evaluate  blast  mitigation  technology.  Historically,  the  Department  of  Army  has  had 
difficulty  collaborating  with  industry  and  academia  on  underbody  blast  events  due  to  the  sensitive 
nature  of  the  work.  Data  generated  from  testing  military  vehicles  is  typically  rated  CLASSIFIED  and  not 
readily  sharable.  To  alleviate  this  issue,  TARDEC  has  fabricated  a  generic  vehicle  hull  (Figure  3.21)  with 
the  intent  to  share  data  with  academia  and  the  industry  to  spur  innovation  in  blast  mitigation 
technologies  [30,  33].The  availability  of  an  existing  LS-DYNA  model  obtained  for  the  Loci-Blast/LS-DYNA 
coupling  work  done  under  SimBRS  WD0034  made  it  an  excellent  choice  for  the  final  tests  of  the  one-way 
coupler.  As  a  first  step  in  the  simulation  process,  an  STL  file  of  the  exterior  of  the  GHULL  LS-DYNA  model 
was  generated  using  LS-PREPOST.  Figure  3.21  shows  that  all  the  details  of  the  external  geometry  of  the 
LS-DYNA  model  are  preserved  in  the  STL  file.  Next,  a  set  of  LOAD_SEGMENT  elements  were  defined  for 
the  undersurface  of  the  GHULL  configuration.  The  parts  of  the  undersurface  that  is  loaded  by  CTH  is 
shown  if  Figure  3.22.  Due  to  the  number  of  elements  in  the  LS-DYNA  model,  this  procedure  resulted  in 
46,729  LOAD_SEGMENT  entries.  This  represents  an  extreme  test  of  the  one-way  coupling  software.  A 
CTH  run  was  created  that  contained  approximately  149  million  cells  (900x360x460  cells).  The  mesh 
spacing  is  a  constant  1cm  in  all  three  axes.  The  soil  was  modeled  using  the  HEP  IMD_INT_SOIL  material 
which  is  an  intermediate  silty  sand  with  a  wet  density  of  2.01  g/cmA3.  The  charge  weight,  depth  of 
burial,  and  position  under  the  vehicle  were  the  same  values  used  in  the  SimBRS  WD0034  simulations. 
CTH  simulations  were  run  out  to  2e-3  seconds  and  the  load  curves  generated  where  included  in  an  LS- 
DYNA  simulation  of  the  complete  vehicle  that  ran  out  to  5e-3  seconds.  For  the  CTH  simulations,  half- 
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symmetry  was  employed  but  the  full  geometry  was  used  in  the  LS-DYNA  runs.  This  tested  a  feature  of 
the  one-way  coupler  that  allows  users  to  input  a  full  LS-DYNA  geometry  into  CTH  but  assume  half  or 
quarter  symmetry  in  the  CTH  simulation.  Results  from  the  simulation  are  given  in  Figures  3.23-3.28. 
Figures  3.23  and  3.24  show  the  displacement  of  the  external  shell  at  5e-3  and  2e-3  seconds.  Figure  3.25 
shows  the  external  pressures  remaining  on  the  underside  of  the  vehicle  at  5e-3  seconds.  Figure  3.26  and 
3.27  show  the  effective  stress  and  displacement  of  the  internal  frame  and  false  floor  at  le-3  and  5e-3 
seconds.  These  results  are  similar  to  those  obtained  in  WD0034.  The  success  of  this  simulation 
demonstrates  that  the  one-way  coupler  can  handle  complex  geometries  and  provide  acceptable 
accuracy  for  analyzing  the  effects  of  buried  mines  on  realistic  vehicle  geometries. 


Figure  3.  21:  TARDEC  Generic  Hull  STL  geometry  from  LS-PREPOST  (above),  physical  hull  (below) 
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Figure  3.  22:  Geometry  with  surfaces  loaded  by  CTH  pressures  shown  in  white 


Figure  3.  23:  Contours  of  external  hull  displacement  at  t=5e-3  seconds 
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Figure  3.  24:  Contours  of  internal  hull  displacement  at  t=2e-3  seconds 


Figure  3.  25:  Contours  of  external  hull  pressure  at  t=5e-3  seconds 
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GHULL  with  CTH  Load  Curves 

Time  =  0.0009995 

Contours  of  Effective  Stress  (v-m) 

max  IP.  value 

min=219372,  at  elem#  604751 
max=3.17828e+08,  at  elem#  608805 


Fringe  Levels 

3.178e+08  _ 

2.861e+08 

2.543e+08 

2.225e+08 

1.908e+08 

1.590e+08 

1.273e+08 

9.550e+07 

6.374e+07 

3.198e+07 

2.194e+05 


Figure  3.  26:  Contours  of  effective  stress  on  the  internal  frame  at  t=le-3  seconds 


Figure  3.  27:  Contours  of  displacement  on  the  internal  frame  at  t=5e-3  seconds 
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Section  4:  Summary  and  Conclusions 

•  Under  SimBRS  WD0037,  two  enhancements  to  CTH  (Version  9.1)  were  made  to  incorporate  the 
ERDC  HEP  geo-material  model  and  provide  a  one-way  coupling  procedure  for  linking  CTH 
generated  pressures  to  an  LS-DYNA  structural  dynamics  simulation. 

•  The  CTH  HEP  model  was  validated  for  buried  explosive  cases  using  2D  witness  plate  simulations 
that  compared  results  obtained  by  CTH  HEP  with  results  obtained  by  Kerley  and  Littlefield  using 
CTH  and  EPIC. 

•  The  CTH  HEP  model  implementation  was  also  used  to  simulate  the  ERDC  IMD  plate  experiments. 
Correlations  with  experimental  data  for  the  three  IMD  soils  indicated  that  CTH  HEP  provides 
excellent  results  for  the  wetter  materials. 

•  Results  for  the  IMD  DRY  SAND  material  illustrated  the  effect  of  the  strength  model  on  impulse 
and  crater  geometry  and  pointed  to  the  need  for  a  more  effective  strength  model  for  the  CTH 
HEP  implementation. 

•  A  detailed  model  of  the  IMD  test  configuration  was  also  run  for  the  WET  CLAY  material  and 
compared  with  a  fit  based  on  the  procedure  developed  by  Kerley.  The  impulse  from  these 
simulations  was  in  excellent  agreement  with  the  experiment. 

•  All  of  the  tests  and  correlations  performed  to  date  indicate  that  the  current  implementation  of 
the  HEP  model  in  CTH  provides  an  accurate  alternative  to  existing  models  for  most  soil  types. 
However,  more  work  needs  to  be  done  to  improve  the  prediction  of  the  effects  of  material 
strength  for  dry  materials. 

•  The  one-way  coupler  was  implemented  and  verified  using  small  witness  plate  simulations  as 
well  as  large-scale  simulations  of  the  DRDC  plate  experiment  and  the  Generic  Hull  experiment. 
The  results  of  these  simulations  indicate  that  the  one-way  coupler  has  the  potential  to  be  an 
effective  tool  for  analysis  of  blast  effects  of  buried  mines  on  vehicle  underbodies. 


Section  5:  Recommendations 

Both  the  HEP  geo-material  soil  model  and  the  one-way  coupler  to  LS-DYNA  have  significant  promise  in 
the  field  of  Underbody  Blast  modeling  for  lEDs.  Evaluations  should  be  continued  to  assess  and  improve 
these  features,  as  well  as  incorporate  them  into  future  versions  of  CTH,  including  the  production  version 
released  by  Sandia  National  Labs.  It  is  recommended  that  the  the  Appendices  in  this  report  be  used  as 
the  starting  point  for  the  routines  that  need  to  be  modified. 

5.1:  Technology  Transition  Opportunities  &  Drops 

These  enhancements  to  CTH  v9.1  are  available  to  any  DoD  agency  or  industry  partner  with  the  following 
restriction.  The  CTH  software  is  owned  and  maintained  by  Sandia  National  Labs,  NM  and  has  been 
designated  Export  Control  Information  (ECI)  under  the  International  Traffic  in  Arms  Regulations  (ITAR). 
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5.2:  Payoffs 

These  enhancements  to  CTH  v9.1  are  beneficial  for  the  following  reasons: 

•  One-way  coupler  with  LS-DYNA: 

o  CTH  has  a  vast  library  of  charge  and  soil  models,  and  can  be  used  to  describe  the  blast 
loading  into  the  structure,  assuming  a  rigid  interface  at  the  underbody. 

o  Once  the  blast  loading  is  determined,  it  is  applied  to  the  Lagrangian  vehicle,  which  runs 
much  faster  than  the  Arbitrary  Langangian-Eulerian  (ALE)  formulation  involving  fluid- 
structure  interaction. 

o  One-way  couplers  have  been  shown  to  be  reasonably  accurate  when  coupled  to  the  full 
two-way  simulations,  considering  the  reduced  computational  costs  for  the  one-way 
treatment. 

•  Development  of  Hybrid-Elastic-Plastic  (HEP)  soil  model: 

o  The  HEP  Geo-material  model  from  ERDC  has  proven  to  provide  accurate  simulations  of 
the  response  of  a  variety  of  geo-materials  to  extreme  loading  due  to  blast,  hypervelocity 
penetration,  and  ground  shock. 

o  Accurate  material  models  for  soils  and  other  geo-materials  are  critical  for  analysis  of 
shallow  buried  explosives  such  as  mines  and  improvised  explosive  devices. 

o  A  large  library  of  HEP  models  has  been  developed  by  ERDC  and  can  be  leveraged  for  use 
in  CTH  directly  and  also  in  the  one-way  coupler. 
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Section  7:  Disclaimer 

Reference  herein  to  any  specific  commercial  company,  product,  process,  or  service  by  trade  name, 
trademark,  manufacturer,  or  otherwise  does  not  necessarily  constitute  or  imply  its  endorsement, 
recommendation,  or  favoring  by  the  United  States  Government  or  the  Dept,  of  the  Army  (DoA).  The 
opinions  of  the  authors  expressed  herein  do  not  necessarily  state  or  reflect  those  of  the  United  States 
Government  or  the  DoD,  and  shall  not  be  used  for  advertising  or  product  endorsement  purposes. 


UNCLASSIFIED:  Distribution  Statement  A.  Approved  for  Public  Release 


68  |  P  a  g  e 


ENHANCEMENTS  AND  ANALYSIS  OF  CTH  SOFTWARE  FOR  UNDERBODY  BLAST 


Section  8:  References/Bibliography 

1.  McGlaun,  J.,  Thompson,  S.  and  Elrick,  M.  "CTH:  A  Three-Dimensional  Shock-Wave  Physics  Code.", 
International  Journal  of  Impact  Engineering,  vol.  10,  pp.  351-360, 1990. 

2.  Zimmerman,  H.,  Wagner,  M.,  Carney,  J.,  Ito,  Y.,  "Effects  of  Site  Geology  on  Ground  Shock 
Environments,  Report  1,  Constitutive  Models  for  Materials  12,  13,  and  W1-W10,"  Technical  Report  SL-87- 
19,  U.S.  Army  Engineer  Waterways  Experiment  Station,  Vicksburg,  MS  1987 

3.  Akers,  S.  and  Adley,  M.,  "Constitutive  Models  Used  to  Simulate  Penetration  and  Perforation  of 
Concrete  Targets,"  Proceedings  of  ASME  Pressure  Vessels  and  Piping  Conference:  Structures  under 
Extreme  Loading  Conditions,  PVP-vol.325,  Montreal,  Canada,  1996. 

4.  Akers,  S.,  Ehrgott,  J.,  Rickman,  D.,  Williams,  E.,  Danielson,  K.,  "Blast  Environments  From  Surface-Laid 
and  Buried  Bare  Charges:  A  Joint  Experimental  and  Simulation  Effort" 

5.  Moral,  R.,  Danielson,  K.,  and  Ehrgott,  J.,  "Tactical  Wheeled  Vehicle  Survivability:  Comparison  of 
Explosive-Soil-Air-Structure  Simulations  to  Experiments  Using  the  Impulse  Measurement  Device,"  Army 
Corp  of  Engineers  Technical  Report  ERDC/GSL  TR-10-27,  August,  2010. 

6.  Danielson,  K.,  Ehrgott,  J.,  Moral,  R.,  Adley,  M.,  O'Daniel,  J.,  Gerlach,  C.,  Johnson,  G.,  Littlefield,  D., 
"Lagrangian  Meshfree  Methodologies  for  Predicting  Loadings  from  Buried  Munitions  Detonations," 
Paper  FO-002,  27th  Army  Science  Conference,  Orlando,  FL,  Nov.  29-Dec.2, 

7.  Tillotson,  J.,  "Metallic  Equations  of  State  for  Hypervelocity  Impact,"  Technical  Report  GA-3216, 
General  Atomics,  San  Diego,  CA  1962. 

8.  Hertel,E.  and  Kerley,  G.,  "CTH  Reference  Manual:  The  Equation  of  State  Package,"  SAND98-0947, 
Sandia  National  Laboratory,  Albuquerque,  NM  1998 

9.  Hertel,  E.  and  Kerley,  G,  "CTH  EOS  Package:  Model  Installation  Guide,"  SAND98-0946,  Sandia  National 
Laboratory,  Albuquerque,  NM  1998 

10.  Brannon,  R.  and  Wong,  M.,"MIG  0.0  Model  Interface  Guidelines:  Rules  to  Accelerate  Installation  of 
Numerical  Models  Into  any  Compliant  Parent  Code,"  SANDIA96-2000,  Sandia  National  Laboratory, 
Albuquerque,  NM  1996 

11.  Akers,  S.  A.,  Private  Communication,  December  2010. 

12.  Kerley,  G.,  "Numerical  Modeling  of  Buried  Mine  Explosions",  ARL-CR-461,  Army  Research  Laboratory, 
Aberdeen  Proving  Ground,  MD  2001 


UNCLASSIFIED:  Distribution  Statement  A.  Approved  for  Public  Release 


69  |  P  a  g  e 


ENHANCEMENTS  AND  ANALYSIS  OF  CTH  SOFTWARE  FOR  UNDERBODY  BLAST 


13.  Kerley,  G.,  "The  Effects  of  Soil  Type  on  Numerical  Simulations  of  Buried  Mine  Explosions/'  Report 
KTS02-3,  Kerley  Technical  Services,  Appomattox,  VA  2002 

14.  Kerley,  G.,  "  On  the  Numerical  Simulation  of  Buried  Mine  Explosions:  Choosing  Constitutive 
Models,",  Report  KTS05-3,  Kerley  Technical  Services,  Appomattox,  VA  2005 

15.  Littlefield,  D.,  "Validation  of  the  Kerley  Soil  Model  in  CTH,"  presented  at  The  81st  Shock  and 
Vibration  Symposium,"  Orlando  FL,  October  24-28,  2010 

16.  De  Vries,  D.,  "Thermal  Properties  of  Soils,"  Chapter  7,  Physics  of  Plant  Environment,  W.R  Van  Wijk, 
Ed.,  North-Holland  Publishing  Company,  Amsterdam,  1966 

17.  Abu-Hamdeh,  N.,  "Thermal  Properties  of  Soils  as  Affected  by  Density  and  Water  Content," 
Biosystems  Engineering,  Vol.  86,  No.  1,  pp.  97-192,  2003 

18.  Hendrickx,  J.,  van  Dam,  R.,  Borchers,  B.,  Curtis,  J.,  Lensen,  H.,  Harmon,  R.,  "Worldwide  distribution  of 
soil  dialectric  and  thermal  properties,"  Proceeding  SPIE  5089,  1158(2003),  Detection  and  Remediation 
Technologies  for  Mine  and  Mine-like  Targets,  VIII. 

19.  Scott,  H.,  Soil  Physics,  Agricultural  and  Environmental  Applications,  Iowa  State  University  Press, 
Ames,  IA,  2002 

20.  Ehrgott,  J.,  "Tactical  Wheeled  Vehicle  Survivability:  Results  of  Experiments  to  Quantify  Aboveground 
Impulse,"  Army  Corp  of  Engineers  Technical  Report  ERDC/GSL  TR-10-7,  March  2010. 

21.  Erhgott,  J,  et  al.,  "Tactical  Wheeled  Vehicle  Survivability:  Construction  and  Characterization  of  Dry 
Sand  Testbeds  in  Bare  Charge  Experiments,"  Army  Corp  of  Engineers  Technical  Report  ERDC/GSL  TR-10- 
15,  March  2010. 

22.  Erhgott,  J,  et  al.,  "Tactical  Wheeled  Vehicle  Survivability:  Construction  and  Characterization  of 
Intermediate-Air-Voids,  Silty  Sand  Testbeds  in  Bare  Charge  Experiments,"  Army  Corp  of  Engineers 
Technical  Report  ERDC/GSL  TR-10-16,  April  2010. 

23.  Erhgott,  J,  et  al.,  "Tactical  Wheeled  Vehicle  Survivability:  Construction  and  Characterization  of  Wet 
Clay  Testbeds  in  Bare  Charge  Experiments,"  Army  Corp  of  Engineers  Technical  Report  ERDC/GSL  TR-10- 
17,  April  2010. 

24.  Williams,  K.,  McClennan,  S.,  Durocher,  R,  St-Jean,  B.,  Trembelay,  J.,  "Validation  of  a  Loading  Model 
for  Simulating  Blast  Mine  Effects  on  Armoured  Vehicles,"  7th  International  LS-DYNA  Users  Conference, 
Detroit,  Ml  2002. 

25.  Bessette,  Greg,  Private  Communication,  November,  2011. 

26.  Silling,  S.A.,  "CTH  Reference  Manual:  Johnson-Holmquist  Ceramic  Model,"  Report  SAND92-0576, 
June,  1996. 


UNCLASSIFIED:  Distribution  Statement  A.  Approved  for  Public  Release 


70  |  P  a  g  e 


ENHANCEMENTS  AND  ANALYSIS  OF  CTH  SOFTWARE  FOR  UNDERBODY  BLAST 


27.  Littlefield,  D.L,  Harvey,  J.P,  and  Kulatha,  S.,  "Implementation  of  a  Uni-directional  Coupling  Interface 
into  CTH,"  HPCMP  PET  Project  Report,  2006. 

28.  Crawford,  D.A.,  "CHSSI  CTH  Presto  One-Way  Coupling  Procedure,"  CTH  Technical  Note,  January, 
2007. 

29.  Thompson,  D.,  Luke,  E.,  Moore,  C.,  Remotique,  M.,  Tong,  X.L.,  Weed,  R.,  and  Ivancic,  P.,  "Enhance 
Simulation  of  Blast-Vehicle  Interactions  using  Loci/BLAST  and  LS-DYNA,",  Final  Report,  SimBRS  EWD  No:, 
WD0034,  June,  2012. 

30.  Dooge,  Dooge  D,  Dwarampudi  R,  Schaffner  G,  Miller  A,  Thyagarajan  R,  Vunnam  M,  Babu  V,  2011, 
"Evolution  of  Occupant  Survivability  Simulation  Framework  Using  FEM-SPH  Coupling,"  2011  NDIA 
Ground  Vehicle  System  Engineering  and  Technology  Symposium  (GVSETS),  MSTV  Mini-Symposium 
Paper,  August  9-11,  Dearborn,  Michigan 

31.  Toussaint,  Genevieve  and  Bouamoul,  Amal,  "Comparison  of  ALE  and  SPH  Methods  for  Simulating 
Mine  Blast  Effects  on  Structures,"  Technical  Report,  DRDC  ValcartierTR  2010-0326,  December,  2010. 

32.  Bessette,  G.  Vaughan,  C.,  Bell,  R.,  Yarrington,  P.,  and  Attaway,  S.,  "Zapotec:  A  Coupled  Eulerian- 
Lagrangian  Computer  Code,  Methodology  and  Users  Manual,  Version  1.0,"  Sandia  Technical  Report, 
Sandia  National  Laboratories,  2003. 

33.  Schaffner,  G  and  Miller,  A,  "Mitigation  of  Occupant  Acceleration  in  a  Mine-Resistant  Ambush- 
Protected  Vehicle  Blast  Event  using  an  Optimized  Dual-Hull  Approach,"  2012  NDIA  Ground  Vehicle 
System  Engineering  and  Technology  Symposium  (GVSETS),  MSTV  Mini-Symposium  Paper,  August  14-16, 
Dearborn,  Michigan. 


UNCLASSIFIED:  Distribution  Statement  A.  Approved  for  Public  Release 


71  |  P  a  g  e 


ENHANCEMENTS  AND  ANALYSIS  OF  CTH  SOFTWARE  FOR  UNDERBODY  BLAST 


Section  9:  Appendices 

APPENDIX  A:  Description  of  HEP  Model  Software  and  CTH  Code 
Modifications 

A.l  Overview 

As  previously  stated  in  this  report,  the  CTH  HEP  EOS  implementation  followed  the  CTH  EOS 
implementation  guidelines  outlined  by  Hertel  and  Kerley  [9].  The  procedure  for  implementing  a  new 
EOS  in  CTH  can  be  summarized  in  the  following  steps: 

1.  Define  a  material  number  for  the  new  material.  The  CTH  HEP  EOS  was  defined  to  be  model  number 
27  in  CTH  9.1 

2.  Modify  CTH  routine  EOSUIP  to  assign  the  number  and  read  the  input  data  specific  to  the  selected 
model  number. 

3.  Create  a  new  routine  to  check  and  verify  the  input  data  and  modify  the  CTH  routine  EOS1SI  to  call  the 
data  check  routine  based  on  the  selected  model  number 

4.  Create  a  new  routine  to  assign  "extra  variables"  needed  by  the  HEP  and  modify  CTH  routine  EOS1SK 
to  call  the  new  extra  variable  routine  based  on  the  selected  model  number 

5.  Create  two  new  driver  routines  for  the  underlying  EOS  routines.  The  first  of  these  driver  routines  is  as 
a  scalar  version  of  the  EOS.  It  takes  in  density  and  temperature  and  returns  pressure,  internal  energy, 
and  sound  speed.  The  scalar  routine  is  called  once  for  each  cell  containing  an  HEP  material.  The  second 
driver  routine  is  a  vectorized  EOS  routine  that  takes  in  density  and  energy  and  returns  pressure, 
temperature  and  sound  speed.  The  vectorized  version  of  the  EOS  driver  takes  in  data  from  an  entire 
array  of  cells  and  processes  the  data  in  a  vector  loop. 

6.  Modify  CTH  routine  EOS  to  call  the  scalar  version  of  the  EOS  driver  and  routine  EOSMRE  to  call  the 
vector  version  of  the  EOS  driver 

7.  Add  a  default  input  block  for  the  new  EOS  to  the  CTH  EOS  data  file.  The  format  for  the  input  block 
requires  a  list  of  keywords  that  will  identify  the  input 

A  similar  procedure  for  implementing  a  constitutive  (strength)  model  in  CTH  is  not  as  well 
documented.  As  with  the  EOS,  the  first  task  is  to  pick  a  material  number.  Unlike  the  EOS,  most  of  the 
modifications  to  existing  CTH  routines  are  implemented  in  the  form  of  several  include  files  in  C  header 
file  format.  These  include  file  perform  a  variety  of  tasks  such  as  setting  dimensions  of  arrays  used  by  the 
input  processor,  defining  extra  variables  and  initial  values,  and  calling  the  main  driver  routine  for  the 
elastic-plastic  strength  model  of  a  particular  material.  The  majority  of  the  include  files  are  used  in  the 
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UINCHK  and  UINEP  routines  that  handle  data  checking  and  data  input.  The  main  elastic-plastic  driver 
routine  for  a  material  is  called  from  ELSG.  The  definition  of  extra  variables  is  performed  in  UINISV. 

The  following  section  lists  both  the  new  routines  that  make  up  the  HEP  EOS  and  strength  model  and 
gives  a  short  description  of  what  they  do  or  the  modification  made  to  an  existing  CTH  routine.  The 
directory  containing  the  source  code  for  the  new  or  modified  routine  is  given  in  parenthesis  after  the 
source  file  name 

A.2  New  EOS  routines 

The  following  routine  are  new  routines  specific  to  the  HEP  EOS  implementation.  Most  of  the  new 
software  described  in  the  following  sections  was  written  in  Fortran  2003  free  field  format  signaled  by 
the  F90  file  extension.  A  new  CTH  subdirectory  (/src/raw)  was  created  to  hold  the  new  HEP  model 
software  along  with  the  "one-way"  coupler  software 

A.2.1  hep_constants_module.F90  (/src/modules) 

Fortran  2003  module  containing  global  data  used  by  the  HEP  model  and  support  routines. 

A.2.2  eshepe.F90  (/src/raw) 

The  "vectorized"  EOS  driver  that  takes  in  density  and  energy  and  returns  pressure,  temperature,  and 
sound  speed.  Called  by  EOSMRE. 

A.2.3  eshept.F90  (/src/raw) 

The  "scalar"  EOS  driver  that  takes  in  density  and  temperature  and  returns  pressure,  energy,  and 
sound  speed.  Called  by  EOS 

A.2.4  eshepi.F90  (/src/raw) 

The  data  check  routine  for  the  HEP  EOS  model.  Called  by  EOS1SI. 

A.2.5  eshepx.F90  (/src/raw) 

The  HEP  EOS  extra  variable  routine.  Called  by  EOS1SK 

A.2.6  hepeos.F90  (/src/raw) 

The  main  HEP  EOS  routine.  Calls  HEPHydroStat  to  compute  energy  independent  component  of 
pressure.  Called  by  eshepe.F90  and  eshepe.F90 

A.2.7  HEPHydroStat.F90 

Computes  the  energy  independent  component  of  the  HEP  pressure.  Called  by  hepeos.F90 
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A.3  CTH  routines  modified  for  HEP  EOS  model 

A.3.1  eosuip.F  (/src/eos) 

Added  CTH  HEP  model  as  model  number  27. 

A.3.2  eoslsi.F  (/src/eos) 

Added  call  to  eshepi.F90  check  routine 
A.3.3  eoslsk.F  (/src/eos) 

Added  call  to  eshepx.F90  extra  variable  routine 
A.3.4  eos.F  (/src/eos) 

Added  call  to  eshept.F  scalar  EOS  driver 
A.3.5  eosmre.F  (/src/eos) 

Added  call  to  eshepe.F  vector  EOS  driver 

A.4  New  HEP  Elastic-Plastic  Model  Routines 

A.4.1  hepchk.F90  (/src/raw) 

Data  check  routine  for  HEP  EP  model.  Called  in  uinchk.F 

A.4.2  hepepx.F90  (/src/raw) 

Extra  variable  routine  for  HEP  EP  model.  Called  by  uinisv.F 

A.4.3  hepepdr.F90  (/src/raw) 

Driver  routine  for  HEP  EP  routines.  Calls  SHEPYieldFunction  and  SHEPShearModulus.  Called  by  elsg.F 
A.4.4  SHEPYieldFunction. F90  (/src/raw) 

Computes  SHEP  yield  strength.  Called  by  hepepdr.F90 
A.4.5  SHEPShearModulus.F90  (/src/raw) 

Computes  SHEP  shear  modulus.  Called  by  hepepdr.F90 

A.5  CTH  routines  modified  for  HEP  EP  model 

A.5.1  uinchk.F  (/src/reb) 

Added  m27pl.h  include  file 
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A.5.2  uinep.F  (/src/reb) 

Added  m27p0a.h,  m27p0b.h,  m27p0c.h,  and  m27p0d.h  include  files. 

A.5.3  uinisv.F  (/src/sas) 

Added  m27p2.h  include  file. 

A.5.4  elsg.F  (/src/sas) 

Added  m27p3.h  include  file 
A.5.5  elygp.F  (/src/sas) 

Added  modifications  to  reference  SHEP  yield  and  shear  modulus  values. 

A.6  New  HEP  EP  model  include  files  in  /include/cth 

A.6.1  m27p0a.h 

Sets  material  model  number,  material  name,  title  and  MIG  routine  array  dimensions. 

A.6.2  m27p0b.h 

Contains  call  to  MIG  RDATA  routine 

A.6.3  m27p0c.h 

Contains  call  to  MIG  LOOKFK  routine 

A.6.4  m27p0d.h 

Initializes  default  EP  routine  values  such  as  TMELT  and  YLDVM 

A.6.5  m27pl.h 

Contains  call  to  hepchk.F  data  check  routine 

A.6.6  m27p2.h 

Contains  call  to  hepepx.F  extra  variable  routine. 

A.6.7  m27p3.h 

Contains  call  to  hepepdr.F  HEP  EP  model  driver  routine 
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APPENDIX  B:  Description  of  CTH-LS-DYNA  coupler  software 
and  CTH  code  modifications 

B.l  Overview 

The  CTH-LS-DYNA  coupler  was  initially  based  on  the  software  developed  by  Littlefield  et  al.  [21]. 
Over  the  course  of  the  development  of  the  current  coupler  much  of  this  software  was  either  highly 
modified  or  completely  rewritten.  However,  there  are  a  few  routines  that  are  have  minor  modifications 
or  are  unchanged  form  the  original  code.  Most  of  the  new  code  for  the  coupler  deals  with  processing 
LS-DYNA  LOAD_SEGMENT  and  LOAD_SHELL  keyword  inputs  in  CTH,  extracting  time  averaged  pressures 
for  a  CTH  rigid  material  and  saving  them  to  a  history  file  named  femh.  Other  routines  handle  reading 
and  writing  restart  files  and  post-processing  the  data  in  the  femh  file  to  generate  a  file  containing  LS- 
DYNA  DEFINE_CURVE  input.  Although  the  new  code  is  focused  solely  on  LS-DYNA  support,  the  original 
code  also  supported  DYNA3D.  Some  of  the  logic  for  DYNA3D  support  was  maintained  but  has  never 
been  tested. 

A  typical  CTH  run  using  the  coupler  software  performs  the  following  tasks 

1.  The  fem/endfem  block  in  the  users  main  CTH  input  file  is  read  by  UINFEM  and  various  flags  and  values 
are  set  based  on  the  users  input. 

2.  The  femh  file  is  created  if  it  doesn't  already  exist. 

3.  The  users  LS-DYNA  input  defining  LOAD_SEGMENT  or  LOAD_SHELL  inputs  is  read  and  used  to  define 
CTH  cell  locations  closest  to  the  target  CTH  LOAD_SEGMENT  or  LOAD_SHELL 

4.  As  CTH  runs  to  the  prescribed  stop  time,  the  coupler  software  computes  a  time  average  pressure  to 
be  applied  to  the  target  LS-DYNA  LOAD_SEGMENT  and  LOAD  _SHELL  element  faces.  The  algorithm  and 
some  of  the  code  that  performs  this  task  was  borrowed  from  the  existing  CTH  rigid  material  model  and 
modified  for  use  by  the  coupler. 

5.  The  computed  pressures  are  saved  in  the  femh  file  at  a  user  defined  frequency. 

6.  At  the  end  of  the  job,  the  femh  file  is  scanned  in  a  post-processing  step  and  the  file  containing  the  LS- 
DYNA  DEFINE_CURVE  input  is  generated. 

Steps  3  and  4  are  performed  in  parallel  for  mpicth  runs.  All  other  steps  are  done  by  the  parent  process. 
The  current  version  of  the  software  is  limited  to  LOAD_SEGMENT  and  LOAD_SHELL  inputs.  This  can  lead 
to  extremely  large  DEFINE_CURVE  files.  An  alternate  approach  that  would  allow  the  use  of 
LOAD_SEGMENT_SETS  to  reduce  the  number  of  generated  load  curves  is  being  investigated  but  is  not  in 
the  current  release  of  the  software.  The  following  sections  describe  the  new  software  and  the  CTH 
routines  that  were  modified  to  support  the  coupler.  The  directory  containing  the  new  software  is  shown 
in  parentheses. 
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B.2  New  CTH-LS-DYNA  coupler  software 

B.2.1  femesr.F/femsew.F  (/src/raw) 

Routines  to  read  and  write  coupler  restart  data.  Called  by  dbmwtf.F 

B.2.2  femop.F  (/src/raw) 

Executive  routine  for  opening  coupler  input  and  output  files 

B.2.3  femini.F  (/src/raw) 

Initialization  routine 
B.2.4  fempst.F  (/src/raw) 

Performs  post  processing  of  femh  file  data  and  generation  of  DEFINE_CURVE  data 

B.2.5  femsur.F90  (/src/raw) 

Computes  surface  normals  and  areas  of  target  LS-DYNA  LOAD_SEGMENT  and  LOAD_SELL  surfaces 
B.2.6  femwrt.F90  (/src/raw) 

Writes  data  to  the  femh  file  at  a  user  defined  frequency 

B.2.7  setfem.F90  (/src/raw) 

Selects  type  of  coupling.  Original  version  of  coupler  had  support  for  DYNA3D.  Kept  for  possible 
future  extension  of  the  code  to  reinstate  DYNA3D  support 

B.2.8  seting.F  (/src/raw) 

Read  in  data  from  DYNA3D  files  to  define  load  curves.  Kept  for  possible  future  extension  to  reinstate 
DYNA3D  support 

B.2.9  setlsd.F  (/src/raw) 

Reads  the  LS-DYNA  input  keyword  file.  Parses  files  for  LOAD_SEGMENT,  LOAD_SHELL  and  NODE 
keywords.  Computes  array  dimensions  and  allocate  memory  for  coupler  arrays. 

B.2.10  settra.F  (/src/raw) 

Finds  location  of  CTH  cell  centers  nearest  to  the  target  LS-DYNA  surface 

B.2.11  hfeme.F  (/src/raw) 

Database  manager  for  hfemeds  routine 

B.2.12  hfemeds.F  (/src/raw) 
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Extracts  pressures  for  rigid  target  material  and  computes  time  averaged  pressure 

B.2.13  uinfem.F  (/src/raw) 

Reads  fem/endfem  block  input  and  sets  flags  needed  by  coupler 

B.2.14  countLoadSegs.F90,  countLoadShells.F90,  countNodes.F90  (/src/raw) 

Routines  that  count  the  number  of  LOAD_SEGMENT,  LOAD_SHELL,  and  NODE  inputs.  These  numbers 
are  used  for  allocating  memory  for  coupler  arrays. 

B.2.15  parseLoadSegs.F90,  parseLoadShells.F90,  parseNodes.F90  (/src/raw) 

Reads  LS-DYNA  LOAD_SEGMENT,  LOAD_SHELL,  and  NODE  inputs 

B.2.16  checkFemSym.F90  (/src/raw) 

Checks  to  see  if  I  and  J  direction  symmetry  flags  where  set  in  fem/endfem  block  and  adjust  node  x 
and  y  values  defined  in  the  LOAD_SEGMENT,  NODE,  and  LOAD_SHELL  input. 

B.2.17  findNodelndex.F90  (/src/raw) 

Provides  a  map  between  node  and  element  indices 

B.2.18  getShelllndex.F90  (/src/raw) 

Gets  index  for  shell  elements 

B.3  CTH  routines  modified  for  the  CTH-LS-DYNA  coupler 

B.3.1  dbmwtf.F  (/src/dbm) 

Calls  femesr  and  femesw  restart  data  I/O  routines 

B.3.2  celop.F  (/src/jmm) 

Calls  fempst  at  the  end  of  a  CTH  run.  Also,  main  CTH  subroutine  that  controls  the  global  solution 
process. 

B.3.3  edende.F  (/src/jmm) 

Calls  femwrt  at  the  end  of  each  step  to  write  data  to  femh 

B.3.4  uins.F  (/src/jmm) 

Calls  uinfem  to  read  fem/endfem  input  block  data.  Also  calls  setfem  to  define  type  of  coupling. 

B.3.5  erdrv.F  (/src/rlb) 
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Calls  femfop,  hfeme  and  femweb 

B.3.6  gcntrl.F  (/src/slt) 

Calls  uinfem  and  setfem.  Part  of  the  old  mpigen  software 
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